AQMB  LBCniM  SMniM 

Int^irated  Design  Analysis 
and  Optimisation  of 
AircnA  Structures 

(LAnalyse  Int^rale  de  la  Conception  et 
roptimisation  des  Structures  des  Aeronefe) 


DTIC. 

ELECTE 


AUGO  5  1992 

A 


I 


This  rmuenal  in  this  pubUaakm  was  assenAkd  to  s>9port  a  Ledure  Series 
under  Ow  sponsorship  of  ^Structum  and  hkaerialsPmd  of  AGARD 
and  the  (^nsukarU  and  ExekanfgtPtopanum  of  AGARD  presented  on 
8di-9thJune iminPasadena,  CA,  UtdeedSmes,  22mh-23nlJune  1992 
in  Lisbon,  Fortugsl  and  2^h—26dtJtme  1992 in  Londtm,  Urdted  Kingdom. 


I 


Oltiribiillen  tnd  Autdnbgbf  on  Bsek  Oasnr 


AGARD-LS-186 


ADVISORY  GROUP  FOR  AEROSPACE  RESEARCH  &  DEVELOPMENT 

7  RUE  ANCELLE  92200  NEUILLY  SUR  SEINE  FRANCE 


OTIC 


AGARD  LECTURE  SERIES  186 

Integrated  Design  Analysis 
and  Optimisation  of 
Aircraft  Structures 

(LAnalyse  Integrale  de  la  Conception  et 
rOptimisation  des  Structures  des  Aeronefs) 


1  Accesion  For  ^ 

NTtS 

CRA&l 

DTIC 

Tab 

L  J 

Ur.aiinOL);;Ced 

u 

J  itii  1 

By  ~1 

D-:t' ioutio  ' /  1 

A 

.  Jlla 

r  C-iX'v.3 

Diit 

1  01 

Ml 

- - 

This  material  in  this  publication  was  assembled  to  support  a  Leeture  Series 
under  the  sponsorship  of  the  Structures  and  Materials  Panel  of  AGARD 
and  the  Consultant  and  Exchange  Programme  of  AGARD  presented  on 
8th— 9th  June  1992  in  Pasadena,  CA,  United  States,  22nd— 23rd  June  1992 
in  Lisbon,  Portugal  and  25th— 26th  June  1992  in  London,  United  Kingdom. 


North  Atlantic  Treaty  Organization 
Organisation  du  Traite  de  I’Atlantique  Nord 


92-21237 


The  Mission  of  AGARD 


According  to  its  Charter,  the  mi.ssion  of  AGARL)  is  to  bring  together  the  leading  personalities  of  the  NATO  nations  in  the  fields 
of  science  and  technology  relating  to  aerospace  for  the  following  purposes: 

—  Recommending  effective  ways  for  the  member  nations  to  use  their  research  and  development  capabilities  for  the 
common  benefit  of  the  NATO  community; 


—  Providing  scientific  and  technical  advice  and  assistance  to  the  Military  Committee  in  the  field  of  aerospace  research  and 
development  (with  particular  regard  to  its  military  application); 

—  Continuously  stimulating  advances  in  the  aerospace  sciences  relevant  to  strengthening  the  common  defence  posture; 

—  Improving  the  co-operation  among  member  nations  in  aerospace  research  and  development; 

— -  Exchange  of  scientific  and  technical  information; 

—  Providing  assistance  to  member  nations  for  the  purpose  of  increasing  their  scientific  and  technical  potential; 

—  Rendering  scientific  and  technical  assistance,  as  requested,  to  other  NATO  bixiies  and  to  member  nations  in  connection 
with  research  and  development  problems  in  the  aerospace  field. 


The  highest  autht)rity  within  AGARD  is  the  National  Delegates  Board  consisting  of  officially  appointed  .senior  representatives 
from  each  member  nation.  The  mission  of  AGARD  is  carried  out  through  the  Panels  which  are  composed  of  experts  appointed 
by  the  National  Delegates,  the  Consultant  and  Exchange  Programme  and  the  Aerospace  Applications  Studies  Programme.  The 
results  of  AGARD  work  are  reported  to  the  member  nations  and  the  NATO  Authorities  through  the  AGARD  series  of 
publications  of  which  this  is  one. 

Participation  in  AGARD  activities  is  by  invitation  only  and  is  normally  limited  to  citizens  of  the  NATO  nations. 


The  content  of  this  publication  has  been  reprinlueed 
directly  from  material  supplied  by  AGARD  or  the  authors. 


Published  May  19V2 

Copyright  ©  AGARD  IM92 
All  Rights  Reserved 

ISBN  92-«.^.‘i-()675-K 


I’riiiti’il  by  SpeauUsed  I'rinlinjf  Sen  ices  Limited 
dtH  'hiftwvll  Lane,  l.ou^hton.  Lsxex 


w 


Abstract 


There  is  a  lack  of  precise  information  on  the  effectiveness  of  specific  methods  in  generating  optimum  designs  for  realistic  aircraft 
structures.  In  this  situation  it  is  difficult  for  designers  to  make  decisions  on  which  systems  to  employ  for  a  given  design  problem 
and  which  developments  to  pursue.  Thus  it  is  neces.sary  for  designers  to  be  aware  of  the  relative  merits  of  the  different  methods 
currently  used  for  the  design  optimi.sation  of  advanced  aircraft. 

This  Lecture  Series  covers  the  methods  available  for  the  computer  based  design  analysis  and  design  optimisation  of  aircraft 
structures.  The  Lecture  Series  deals  with  the  principles  and  practices  adopted  to  integrate  the  various  factors  which  are 
considered  in  the  design  of  advanced  aircraft.  These  factors  include:  structural  shape,  aerodynamics,  active  control  technology 
and  aircraft  performance.  Realistic  case  studies  are  used  to  illustrate  the  methods  u.sed  for  different  design  problems. 

The  following  topics  are  covered  in  detail: 

—  Overview  of  integrated  design  analysis,  background,  methods,  objectives  and  requirements. 

—  Optimi.sation  in  design  (CAE/CAD). 

—  A  system  approach  to  aircraft  optimisation. 

—  Case  studies  for  different  design  problems. 

This  Lecture  Series,  sponsored  by  the  Structures  and  Materials  Panel  of  AGARD,  has  been  implemented  by  the  Consultant  and 
Exchange  Programme. 


Abrege 


11  y  a  un  manque  d'informations  precises  sur  I'efficacite  des  methodes  specifiques  qui  ont  ehe  elaborees  pour  I'optimisation  des 
etudes  en  vue  de  la  realisation  de  structures  d’aeronefs.  Dans  cette  situation  il  est  difficile  pour  les  concepteurs  d'avion  de 
decider  des  systemes  a  employer  pour  resoudre  tel  ou  tel  probleme  de  conception  et  d'identifier  les  developpements 
interessants.  II  importe  done,  de  sensibiliser  les  concepteurs  sur  la  valeur  relative  des  differentes  methrrdes  employees  pour 
I'optimisation  de  la  conception  des  aeronefs. 

Ce  cycle  de  conferences  couvre  les  methodes  disponibles  pour  I'analyse  de  la  conception  assistee  par  ordinateur  et 
I'optimisation  de  la  conception  des  structures  d'aeronefs,  II  examine  les  principes  et  les  pratiques  adoptes  pour  I'integration  des 
diffcTents  facteurs  pris  en  compte  lors  de  la  conception  des  aeronefs.  Ces  facteurs  comprennent:  la  forme  structurelle. 
I'aerodynamique,  la  technologic  des  eommandes  actives  et  les  performances.  Des  etudes  de  cas  reelles  sont  utilisees  pour 
illustrer  les  methodes  employees  pour  resoudre  divers  problemes  de  conception. 

Les  questions  suivantes  sont  examinees  dans  le  detail: 

—  panorama  de  I'analyse  integree  de  la  conception,  historique,  methodes,  objectifs  et  besoins 

—  i'optimisation  de  la  conception  (lAO/CAO) 

—  une  approche  “systemes"  a  I'optimisation  des  aeronefs 

—  des  etudes  de  cas  pour  des  problemes  de  conception. 

Ce  cycle  de  conferences  est  presente  par  le  Panel  AGARD  des  structures  et  materiaux;  et  organise  dans  le  cadre  du  programme 
des  Consultants  et  des  Echanges. 
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1.  Introduction 

Structural  Optimisation  is  concerned  with  the  computerised  automatic  design 
of  structures  which  are  optimum  with  respect  to  some  major  design  parameter. 
In  the  aircraft  industry  this  parameter  has  usually  been  structural  weight, 
though  cost,  performance  or  other  factors  are  now  being  considered.  The 
parameter  being  optimised  is  referred  to  as  the  objective  function  and  the 
variables  which  can  be  changed  to  achieve  the  desired  optimum  are  referred  to 
as  design  variabies.  Mathematically  this  can  be  characterised  by  saying  that 
the  problem  is; 

minimise  (or  maximise)  f(x)  x  6  R" 

subject  to  the  constraints  g^(x)  t  0  j  =  1  .  ni 

h^(x)  =  0  k  =  1  .  p 

where  the  design  variables  x  e  R"  are  positive  and  the  range  of  x  for  which 
the  constraints  are  not  violated  constitute  the  feasible  region.  If  the 
objective  function  f(x)  is  structural  weight  the  design  variables  are  size 
parameters  such  as  bar  cross-sections,  plate  thicknesses  and,  in  certain 
cases,  shape  parameters  which  vary  the  geometrical  configuration  of  the 
structure.  Current  researches  are  seeking  to  extend  the  scope  of  structural 
optimisation  to  cover  more  extensive  objective  functions  which  include 
factors  such  as  performance,  cost,  etc.  Indeed,  certain  commercially 
available  systems  already  cover  non-weight  objective  functions.  The 
constraints  on  the  optimum  will  include  behavioural  parameters  so  that  the 
terms  g  (x)  could  include,  stress,  displacement,  flutter  speed,  vibration 

limits  or  any  other  relevant  parameters.  In  addition,  to  behavioural  aspects 
these  constraints  aiso  cover  physical  limits  imposed  by  practical 
manufacturing  considerations  such  as  gauge  limits.  Whilst  equality 
constraints  are  not  common  in  minimum  weight  optimisation  they  can  occur 
where  design  codes  are  empioyed  or  where  components  can  be  selected  from  a 
specific  remge  (i.e.  stock  items). 
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In  all  this  variety  two  aspects  remain  constant  in  all  current  structural 
optimisation  applications.  First,  the  general  problem  which  is  chemacterised 
by  (1.1)  remains  unchanged  so  that  the  basic  nature  of  the  optimisation 
problem  is  the  same  for  all  applications.  Thus,  the  theory  described  in  this 
Lecture  Series  can  be  used  in  all  design  applications.  Secondly,  the 
structural  behaviour  of  the  optimisation  problem  is  always  characterised  by 
the  Finite  Element  Method.  In  many  cases  this  has  lead  to  the  development  of 
optimisation  modules  which  form  an  integral  part  of  many  commercially 
supported  F.E.  packages  i.e.  NASTRAN,  SAMCEF,  IDEAS,  ELFINI,  ANSYS,  etc.  In 
addition,  independent  structural  optimisation  systems  have  been  developed, 
such  as  the  DRA/SCICON  STARS,  MBS  Lagrange,  systems,  which  can  be  attached  to 
any  existing  FE  system.  These  developments  have  resulted  in  Structural 
Optimisation  Methods  being  routinely  available  to  users  of  modern  CAD 
systems. 

The  use  and  application  of  these  methods  in  a  safe  and  effective  manner 
requires  some  understanding  of  the  underlying  mathematical  principles.  As  in 
the  case  of  the  Finite  Element  Method  the  basic  mathematics  provides  a 
’tool-kit*  which  is  repeatedly  used  to  develop  solution  methods.  It  is  shown 
in  later  sections  that  this  process  of  developing  solution  methods  use  the 
optimisation  criteria  as  the  basis  for  creating  the  up-date  formulae  which 
are  the  solution  algorithm  drivers.  Thus,  this  first  part  of  the  Lecture 
Series,  describes  the  optimality  criteria,  the  associated  duality  theory  and 
the  algorithms  themselves. 


2.  A  Basic  Algorithm 

The  computer  based  numerical  solution  process  for  the  problem  defined  at 
(1.1)  is  in  essence,  simple.  It  requires  that  a  repetitive  formu^^e  is^  used 
which  starts  with  an  initial  estimate  of  design  variables  x'”  e  R  and 

systematically  changes  them  until  a  set  x'^  is  generated  after  I  iterations 
which  satisfy  (1.1). 

The  process  is  best  demonstrated  by  considering  an  optimisation  problem  which 
has  no  constraints,  thus  we  seek  to  solve  the  problem 

minimise  f(x)  x  e  R" 

The  solution  can  be  found  using  the  following  solution  algorithm: 


Basic  Algorithm. 

Step  1 

Select  starting  values 

Step  2 

Set  k  =  0 

Step  3 

Set  k  =  k  +  1 

Step  4 

Set  x**'*"  -  Atx'"’) 

f“-'’  =  f(x“"‘’) 

and  choose  a  value  for  € 
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Step  5  If  B(x*'‘***)  2  €  go  to  Step  6; 

Otherwise  go  to  Step  3 

Step  6  Set  x  *  x*"****,  f  =  f(x*'‘***); 

STOP. 

Thus,  step  4  generates  a  new  version  of  the  desi^  varjaj&les  x*”'***  from  the 
earlier  values  x'“  and  the  formula  for  doing  this  A(x‘  J  is  k^^gwn  as^the 
up-date  formula.  I^ause  we  are  moving  from  one  position  x‘“  in  R"  to 
another  position  x'“*  in  r"  AIx**’)  this  constitutes  a  move  of  specific 
length  along  a  given  direction.  Hence 

A(x“‘’)  =  x‘«  ^  a  M  =  x"'*” 

where  p  is  a  direction,  in  r",  from  x****  and  a  is  a  value  giving  the  distance 

to  be  moved  in  this  direction.  Two  questions  now  need  to  be  answered: 
what  should  be  used  for  p 

how  far  to  move  along  p  -i.e.  what  is  the  value  for  a? 

In  answering  the  first  of  these  questions  we  must  select  a  direction  which 
points  towards  the  optimum  and  one,  very  effective,  method  for  achieving  this 
is  to  enforce  the  satisfaction  of  the  optimality  criterion.  Once  a  direction 
has  been  selected  the  value  of  a  is  found  by  seeking  the  minimum  value  of 
f(x)  along  the  direction  . 

EXAMPLE  (Newton’s  Method) 

Suppose  that  the  function  to  be  minimised  f(x)  has  first  zmd  second 
derivatives  available  so  tjgt  it^can  be  approximated  by  a  second  order  Taylor 
expansion  about  a  point  x'  €  r 

f(x“‘’  +  5x1  =  fCx*"’)  +  V /^(x"’).  5x 

^1-2) 

♦  -  «x\  H(x‘'‘’)6x 
2  -  -  - 

where  3x  is  an  increment  x******  -  x**'*;  Af(x*'‘*),  H(x**'*)  are  the  first 
derivative  and  the  Hessian  respectively  of  f*"^*  at  x***^  €R‘*. 

In  order  to  generate  the  up-date  formula  we  note  that  the  optimising 
condition  for  f(x)  X€R"  is: 

Vf(x)  »  0 

If  this  is  enforced  on  (1.2)  with  6x  as  the  free  variables  (since  the  start 
point  x  is  fixed)  then 
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7  f(x)  =  V  f(x“‘’)  +  =  0 

3x  =  -  ^■‘(x"'’).^  fCx"*) 
and  thus  we  choose  fi  =  -H’^Cx****).?  fCx***'). 

If  f(x)  is  of  higher  order  than  a  quadratic  in  xeR"  then  a  is  found  by 
minimising 

f(x'“*+  ofi)  along  y 


(k) 

Having  discussed  the  generation  of  the  up-date  formula  A(x  )  now  return 
to  the  second  unexplained  term  in  the  algorithm,  namely  Bix**  i.  This  is 
simply  a  stopping  criteria!  Because  the  algorithm  is  a  computational  process 
the  optimum  point  is  oniy  located  to  a  specified  level  of  accuracy.  The  term 
B  is,  therefore,  an  accuracy  measure  and  can  be  represented  by  the  change  in 
objective  function  during  an  iteration,  or  the  design  variables.  As  the  next 
section  shows  a  very  effective  measure  is  given  by  noting  the  difference 
between  the  feasible  value  of  the  objective  function  and  the  associated  dual. 


Although  the  algorithm  described  above  is  simple  in  concept  it  is  applicable 
to  all  optimisation  problems,  the  difference  in  allying  it  to  a  range  of 
problems  lies  in  the  changes  associated  with  A(x‘”).  In  the  next  section 
the  optimality  criteria  for  (1.1)  is  introduced  and  algorithms  developed  from 
it. 


3.  Optimality  Criteria  and  Duality 

In  order  to  generate  the  up-date  formula  in  the  previous  section  the 

optimality  criteria  was  used  to  generate  the  solution  algorithm.  However, 

this  applied  only  for  the  case  of  an  unconstrained  optimisation  problem  which 
does  not  represent  the  situation  described  in  (1.1).  Generating  the 

optimality  conditions  for  the  constrained  problem  is  done  by  sequentially 

moving  from  an  unconstrained  problem  to  an  equality  constrained  one  then, 
f’'  illy,  to  the  full  inequality  constrained  optimisation  which  is  the  heart 
of  the  structural  optimisation  design  problem.  This  is  not  a  complex  process 
but  is  too  lengthy  for  inclusion  here  and  is  fully  described  in  reference 
(11. 


Because  the  inequality  constrained  problem  is  the  most  common  form  for 
structural  optimisation  we  shall  consider  a  reduced  form  for  (1.1)  and,  for 
simplicity  take  the  problem  to  be  a  minimisation.  Thus  (1.1)  becomes: 


minimise  f(x) 

subject  to  g^(x)£0  j  =  1  ....m 


(1.3) 


and  the  optimality  criteria  for  ^this  problem  are  known  as  the  Kuhn-Tucker 
conditions.  These  state  that  x  €R"  is  a  local  optimum  for  (1.3)  if  these 
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exist  XeE™  such  that: 


V  Hx")  -  g(x*)  =  0 

~  /  • 

A^g(x  )  =  0 

XiO  g(x*)fcO 

The  first  part  of  (1.5)  are  the  constrained  derivatives  of  the  objective 
function.  That  is,  the  gradient  of  the  objective  function  projected  onto  the 
linearised  form  of  the  constraints.  In  essence  this  projection  returns  us  to 
an  unconstrained  optimisation  problem  so  that  the  algorithm  developed  in 
section  2  once  more  applies. 

An  alternative  form  for  (1.3)  uses  the  Lagrangian  which  is  defined  by 

L(x,X)  =  f(x)  -  x’^.g(x) 

in  which  case  (1.4)  can  be  re-written  as 

7  L(x*.X)  =  0 

X 


X^.g(x*)  =  0 


(1.5) 


XaO  g(x*)sO 


The  standard  problem  defined  by  (1.3)  is  clearly  a  minimisation  problem  which 
is  often  called  the  ‘primal  problem’.  Associated  with  this  is  maximisation 
problem  known  as  the  'dual  problem’  where  a  new  function  is  maximised  subject 
to  a  new  set  of  constraints.  These  two  problems  are  connected  by  a  saddle 
point  so  that  the  minimum  value  which  represents  the  solution  of  the  primal 
problem  is  also  the  value  which  is  the  maximum  value  for  the  dual  problem. 
The  dual  has,  therefore,  two  uses  both  of  which  have  been  exploited  by  the 
developers  of  structural  optimisation  prograunmes.  First,  the  dual 
formulation  provides  an  alternative  description  of  the  optimisation  problem 
which  can  be  used  to  create  solution  algorithms.  Secondly,  it  has  been 
employed  as  a  method  for  generating  bounds  on  the  optimum  which  can  play  the 
rdle  of  the  accuracy  check  B(x'“'). 

Many  forms  for  the  dual  can  be  developed  and  are  discussed  elsewhere, 
ref.  12).  The  one  usually  employed  for  structural  optimisation  states  that 
the  dual  associated  with  the  primal  problem  (1.3)  requires  that,  for  f(x) 
convex  and  g(x)  concave,  we 


where  7x 


-( 
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maximise 

L{x,X) 

subject  to 

7  L(x,X)  *  0 

X 

(1.6) 

XSO 

Although  not  obvious  from  this  formulation,  for  many  structural  problems,  the 
dual  constraints  can  be  solved  to  yield  x(A)  giving  rise  to  an  unconstrained 
optimisation: 

maximise  L(x(A),A) 

with  x(X)  the  solution  of  V  L(x,X)  =  0.  As  we  shall  see  later  this  allows 

X 

the  creation  of  powerfull  dual  solution  algorithms.  These  have  formed  the 
basis  of  the  class  of  dual  algorithms  successfully  employed  in  the  SAMCEF 
system  for  several  years. 


4.  Structural  Optimisation  Algorithms 
4.1  STOESS  RATIOING 

The  Kuhn-Tucker  conditions  (1.4)  or  (1.5)  provide  the  optimality  criterion 
which  are  used  to  generate  the  up-date  formulae  employed  in  the  modern 

automated  design  systems.  But  a  very  simple  up-date  formula  has  been 

effectively  employed  both  as  a  hand  calculation  and  as  a  computerised 

optimisation  method.  This  assumes  that  the  optimum  is  a  vertex  solution  in 
constrained  design  space.  It  required  that  the  number  of  constraints  in 
(1.3)  is  equal  to  or  exceeds  the  number  of  design  variables  i.e.  m^n.  In 

this  situation  the  optimum  is  found  by  solving  n  equality  constraint 
equations 


g  (x)  =  0  k  =  1  ....n 
k 

and,  used  iteratively,  this  produces  the  classical  stress  ratioing  algorithm. 

For  the  minimum  weight  design  of  statically  determinate  structures  subject  to 
stress  constraints  only,  this  method  will  find  the  optimum  in  a  single 
iteration.  For  indeterminate  structures  there  is  no  guarantee  that  an 
up-date  formula  based  on  enforcing  vertex  solution  will  locate  the  optimum 
design.  This  is  because  the  solution  process  takes  no  account  of  the  desire 
to  minimise  the  structural  weight.  For  problems  involving  constraints  other 
than  stress  the  approach  is  highly  inappropriate.  However,  it  is  robust  and 
does  not  require  the  calculation  of  any  derivatives  so  is  effective  in  the 
initial  stages  of  any  solution  involving  a  structural  optimisation  problem  in 
which  stress  constraints  play  a  rdle.  For  this  reason  the  stress  ratio 
algorithm  is  available  in  all  systems  used  for  the  design  of  minimum  weight 
aircraft  structures. 


4.2  OPTIMALITY  CRITERION  ALGORITHMS 


In  the  past  the  optimality  criterion  method  was  a  term  applied  to  a  set  of 
algorithms  devised  by  Venkayya,  Khott  and  Berke  13]  in  the  United  States  and 
Kerr  [4]  in  the  U.K.  However,  these  methods  are  a  special  case  of  the 

general  process  of  using  the  Kuhn-Tucker  conditions  as  the  basis  of  the 
up-date  formulae.  The  term  Optimality  Criterion  methods  is,  more  properly, 
applied  to  a  general  class  of  solution  methods.  In  this  situation  the 

unifying  factor  is  that  the  constraints  are  always  approximated  by  a  first 
order  Taylor  expansion.  (Though  more  recent  work  has  attempted  to  employ 
second  order  approximations  for  the  constraints,  but  this  is  not  considered 
in  the  present  paper.)  The  differentiating  factor  between  the  various 
methods  in  this  class  is  the  order  of  the  Taylor  expansion  used  to 

approximate  the  objective  function.  A  conunon  factor  to  all  the  methods  is 
the  need  to  select  from  the  total  number  of  constraints  a  sub-set  which  are 
considered  to  be  active.  This  means  that,  at  each  iteration,  the  solution 
process  must  establish  which  constraints  are  good  candidates  for  being  strict 
equality,  as  opposed  to  inequality  constraints,  at  the  optimum  point. 

4.2.1.  Linear  Approximation.  The  first  optimality  criterion  approach  assumes 
that  both  the  objective  function  and  the  constraints  are  approximated  by 
first  order  Taylor  expansions.  Thus,  (1.3)  now  becomes 

minimise  f(x)  +  7  f.6x 

subject  to  g(x)  +  N.6x  a  0 

where  N  is  the  matrix  of  constraint  gradients  7g(x)  taken  with  respect  to  the 
design  variables.  These  gradients  can  be  computed  in  a  variety  of  ways  but 
the  finite  element  method  lends  itself  to  analytic  derivatives  for  a  range  of 
element  types.  For  complex  problems  recourse  may  be  made  to  semi-analytic 
derivatives  or,  if  absolutely  necessary,  to  finite  difference  schemes.  The 
generation  of  gradient  derivatives  is  not  discussed  here  as  it  is  a  well 
documented  procedure  available  in  standard  texts. 

The  Lagrangian  associated  with  the  linearised  problem  is: 

L(x,A)  =  f(x)  +  ~  A\(g(x)  +  N.6x) 

The  differential  part  of  the  Kuhn-Tucker  condition  are  then 

7  L(x.A)  =  7  f  -  n\a  =  0 
The  lagrangian  multipliers  can  now  be  extracted: 

A  =  (N.N’^)'‘.N.7_f 

and  can  be  used  to  generate  the  constrained  derivative. 

•V 

The  up-date  formula  can  now  be  constructed  on  the  basis  that  the  optimum 
can  be  located  along  the  direction  of  steepest  descent.  Thus,  the  algorithm 


7_fc  =  |i  -  (N.n’^)'‘.N 
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uses  (-Vf  )  as  the  direction  u  and  a  line  search  can  then  be  conducted, 
c  max  ~ 

If  a  single  constraint  (displacement)  only  is  active,  the  equations  transform 
into  those  used  by  Knott/Venkayya/Berke  and  Kerr  for  the  original  optimality 
criterion  method.  This  formulation  can  also  be  used  to  generate  a  solution 
algorithm  based  on  the  premise  that  the  optimum  can  be  found  from  a  sequence 
of  linear  programmes. 

4.2.1.  Quadratic/Linear  Application.  The  next  level  in  the  hierarchy  of 
solution  methods  assumes  that  the  objective  function  is  approximated  by  a 
second  order  Taylor  expansion  and  the  constraints,  as  at  4.2.1,  by  a  first 
order  expansion.  Thus,  the  problem  (1.3)  becomes: 

minimise  f(x)  +  7  f^.5  x  +  -  6  x.H.6x 
~  ~  2~~~ 

subject  to  g  (x)  +  N.SxtO 

where  H  is  the  Hessian  of  the  objective  function.  The  Lagrangian  associated 
with  this  problem  is: 

L(x.A)  =  f(x)  +  7_f^.3x  +  £«x.H.ax  -  x’'(g(x)  +  N6x 
Thus,  for  optimality: 


y  Ux.A)  =  y  f  +  H.5x  -  n\a  =  0 

giving 

5x  =  -H'‘.(y  f  -  N^A) 

Noting  that  the  constraints  involved  in  this  formula  are  the  active  sub-set 
then: 


g^(x)  +  N.5x  =  0 
or 

N.dx  =  -g__(x) 

=  -N.H'*.(y_f  -  N^.A) 


Thus 

A  =  (N.H'‘.N^)'‘.(N.H‘‘.y  f  -  g(x)) 


and  substituting  back  into  the  expression  for  5x  gives 


or 


6x  =  -H‘‘.|^y  f  -  N^.|(N.H''.N’^)'‘.(N.H'‘.y  f  -  g(x))Jj 
5x  =  -H'\n’^.(n\h"‘.N^)' g(x)  ■  |i  "  y’^-^N.H'*.N’^)'‘.N.H'‘J.H'‘.y  f 


As  in  section  2  this  is  a  Newton  step  and  gives  the  direction  p. 


It  has 
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two  components,  the  first 


-  H'\N‘^.(N.H'\N^r^g(x) 

steps  onto  the  active  constraints,  and  the  second 

projects  the  Newton  direction  onto  the  plane  (linearised)  of  the  active 
constraints. 

The  Quadratic/Linear  approximation  is  a  very  popular  algorithm  and  has 
found  application  in  a  variety  of  systems  including  STARS,  OPTISEN,  OPTI, 
ASTROS.  It  is  used  with  direct  design  variables,  inverse  variables  or 
asymptotic  variables.  Also  the  precise  implementation  of  the  algorithm  may 
vary  from  system  to  system  and  a  range  of  generalised  quadratic  programming 
methods  have  been  employed. 

4.2.1.  Active  Set  Strategies.  As  indicated  above  the  algorithms  described 
rely  on  the  fact  that  the  set  of  constraints  being  used  at  each  iteration  of 
the  algorithm  are  a  sub-set  of  the  total.  These  control  the  space  available 
for  the  up-date  search  direction  to  sweep.  In  principle  this  reduced  set  of 
constraints  generate  a  local  feasible  direction.  The  projection  vectors 
developed  in  4.2.2  restrict  the  search  direction  to  lie  along  this  reduced 
set  which  are  actively  controlling  the  up-date  process.  These  constraints 
are,  therefore,  called  the  active  set  of  constraints  and  the  procedure  for 
deciding  on  which  of  the  total  constraints  are  active,  at  any  iteration,  is 
known  as  the  Active  Set  Strategy. 

There  are  two  parts  to  this  strategy;  one  part  deciding  which  constraints 
are  to  be  included,  the  other  deciding  those  to  be  dropped  from  the  active 
set.  The  process  of  including  a  new  constraint  is  straightforward;  at  each 
iteration  an  analysis  must  be  performed  and  should  any  constraint  be  seen  to 
be  violated  it  must  be  included  in  the  active  set.  Dropping  constraints  is  a 
little  more  complicated.  As  shown  in  the  Kuhn-Tucker  optimality  conditions 
the  Lagrange  multipliers  must  be  positive.  Thus,  any  constraint  at  any 
iteration  which  has  a  negative  Lagrange  multiplier  should  not  be  in  the 
active  set.  The  formulae  in  4.2.2  which  calculate  the  Lagrange  multipliers 
can  be  used  to  identify  these  constraints  with  negative  multipliers.  Such 
constraints  can  then  be  dropped  before  the  next  iteration  is  performed. 

The  exact  process  of  implementation  is  a  little  more  complicated  than  the 
outline  given  above  and  anti  zig-zag  rules  need  to  be  imposed  to  provide  a 
degree  of  smoothness  to  the  operation  of  the  strategy.  Nevertheless,  the 
basic  principles  of  most  active  set  strategies  are  those  given  here. 


5.  Exploiting  the  Dual 

In  this  section  the  power  of  the  dual  formulation  is  demonstrated  both  as  the 
basis  for  a  solution  algorithm  and  as  a  bounding  procedure.  It  is  convenient 
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to  be  a  little  more  specific  in  formulating  the  problem  and,  to  this  end.  a 
linear  weight  function  is  taken  as  the  objective  function.  However,  to 
assist  in  the  process  of  linearising  the  constraints  this  is  transformed  into 
a  non-linear  form  by  the  use  of  inverse  sizing  variables.  The  generic 
structural  minimum  weight  design  problem  becomes 

n  Wj 

minimise  V  - 


inimise  Y  — ~ 

V=1  *1 


where  w^  is  the  specific  weight(mass)  associated  with  the  i  th  design 
variable. 

5.1  DUAL  BOUNDING 

The  dual  problem  associated  with  the  primal  optimisation  problem  (5.1)  is 


ImiMUx.i).?  X  x-b] 

T=i  1  j=i  •’  V  1=1  ■> 


maximise 


w.  m 


subject  to 


Multiply  each  of  the  constraint  equations  by  and  sum  i  =  1  .  n  gives 

n  W  m  n 


-  r  — ^  *Y\YcX=0 

V-.  X.  IV.  U  * 


1=1  1  J=1  1=1  ' 

and  substituting  this  into  the  dual  problem  gives  a  new  dual: 


maximise  )  A  b 
V-.  1  1 


in  w 

subject  to  )  A  c  =  -  i  =  1  .  n 

V_.  1  U 


If  we  explicitly  take  into  account  the  positivity  of  the  design  variables 
this  problem  becomes: 
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minimise  )  X  b 
4_.  i  J 


subject  to 


y  X  b 

m 

y  A  c 
U  J 


*s  — —  I  =  1  .  n 

2 

X 

1 


Thus  (5.3)  is  the  linearised  dual  to  the  primal  problem  (5.1)  and  can  be  used 
to  provide  a  pseudo-dual  bound  on  the  optimum. 

In  order  to  demonstrate  the  procedure  we  assume  that  at  the  end  of  the  kth 
iteration  of  any  of  the  algorithms  described  in  section  4  we  have  the  current 


estimate  of  the  "optimising"  inverse  design  variables  x| 
be  fed  into  (5.3)  to  provide  a  ’dual’  problem 


These  may  now 


minimise 


subject  to 


t  X  b 
^  J  J 

J  =  i 

m 

y  X  c 

J  ‘1 


(x;v 


and  because  the  variables  x^  i  =  1  .  n  are  fixed  the  above  is  a  linear 

programming  problem  which  can  be  solved  to  produce  a  set  of  dual  variables 
A  j  =  1  .  m.  These  can  be  used  to  compute  a  value  of  the  Lagrangian 


(k)  ^(k) 


X  )  which  cam  be  compamed  with 


(1 

weight  =  y 


to  give  a  bound  on  the  optimum.  The  gap  between  these  two  values  is  known  as 
the  duality  gap  and  can  be  used  as  the  function  BIx***  i  used  as  a  stopping 
criteria  in  the  basic  algorithm  of  section  2. 

This  approach  of  using  a  linear  programming  routine  to  solve  a  linearised 
dual  to  provide  a  bound  on  the  optimum  is  used  in  many  structural 
optimisation  codes  including  STARS  and  OPTISEN. 

5.2  A  DUAL  ALGORITHM 


In  order  to  create  a  dual  based  algorithm  we  note  that  the  dual  problem  can 
be  developed  to  remove  any  dependence  on  the  design  variables  x**'  i  = 

1  .  n.  This  is  done  by  explicitly  solving  the  dual  constraints  to  obtain 
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This  expression  for  the  design  variables  can  be  substituted  into  the  dual 
objective  function  so  that  (5.2)  becomes,  simply; 

maximise  L(X) 

which  is  an  unconstrained  maximisation  problem  for  which  the  optimality 
criteria  is 

V_L(X)  =  O 

Assuming  a  2nd  order  Taylor  expansion  for  L(X)  gives  a  new  unconstrained 
maximisation  problem 

maximise  L(X)  ♦  7  L(X)5  X  +  -  5  X  H(X)  5  X 
_  ^  2  ^  _ 


for  which  the  optimum  is  given  by 

V_L(X)  =  A_L(X)  +  HU)  5  X  =  0 
or 

5  X  =  -H’V  L 
(k) 

This  provides  the  up-date  formulae  A(X  )  in  terms  of  the  Lagrange 
multipiiers: 

^  ^(k)  _  H'V  L  (5.4) 


We  note  that: 


VL  has  terms  =  )  c.  x  -  b  =  g  (x) 

OA  t  IJ  1  j  j 


and  H  has  terms  -3^ — = 
ax^  3X^ 

a 

ax^ 

n 

y 

dx 

1 

L 

1=1 

ax 

1 

dX^ 

s 

n 

-L 

c 

t  J 

2w 

3 

X. 

The  up-date  formulae  (5.4)  is  used  in  the  module  OPTI  in  the  SAMCXT  system 
and  was  created  by  Fleuxy,  (Ij.  However,  it  was  used  as  the  basis  of  Kerr’s 
up-date  formulae  within  his  ’optimality-criterion’  code  which  still  forms  the 
core  of  the  B.Ae  Warton  code  ECLYPSE. 
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6.  Conclusion 

The  above  sections  show  that  fundamental  principles  upon  which  the  modern 
structural  optimisation  codes  are  based  are  soimd.  The  theory  with  respect 
to  sizing  variables  is  well  established  and  the  following  chapters  exploit  it 

for  a  wide  range  of  applications.  It  also  forms  a  secure  platform  on  which 
current  research  can  build  to  develop  new  applications.  The  use  of  shape 

variables  represents  one  such  development  and  the  inclusion  of  performance 
characteristics  is  another.  Thus  alternative  design  variables  and  objective 
functions  are  being  introduced  but  the  basis  of  algorithm  development  still 
remains  the  same  as  that  given  above. 
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Summary 

The  structural  optimisation  system  MBB-Lagiange  allows 
the  optimisation  of  homogeneous  isotropic,  orthotropic  or  an¬ 
isotropic  structures  as  well  as  fiber  reinforced  materials.  With 
the  simultaneous  consideration  of  different  requirements  in  the 
design  of  aircraft  structures  it  is  possible  to  reduce  the  numoer 
of  iteration  steps  between  design,  analysis  and  manufacturing. 

Based  on  finite  element  methods  for  structures  and  panel 
methods  for  aerodynamics,  the  analysis  with  sensitivity  in¬ 
cludes  modules  for  static,  buckling,  dynamic,  static  aeroelasdc 
and  flutter  calculations. 

The  opdmtsation  algorithms  consists  of  mathematical  pro¬ 
gramming  methods  and  an  optimality  criteria  procedure. 

The  important  link  between  optimisation  and  analy¬ 
sis/sensitivity  is  the  optimisation  model  which  leads  to  a  very 
modular  architecture. 

Typical  application  examples  show  the  power  and  gener¬ 
ality  of  the  approach. 


1.  INTRODUCTION 

Modem  aircrafts  are  complex  systems  whose  performance 
depends  on  the  interaction  of  many  different  disciplines  and 
parts.  The  complexity  ct  the  problems,  that  means  the  cou¬ 
pling  among  a  very  large  set  of  governing  equations,  is  dealed 
traditionally  by  solving  only  a  subset  of  the  system,  such  as 
aerodynamics,  structures,  flightmechanics,  ccmtrols,  etc. 
(Fig.  1.1  [1]).  For  these  individual  disciplines  great  advances 
were  made  due  to  thetnetical,  computational  and  methodol¬ 
ogy  break  throughs.  Howevo',  these  more  sophisticted  meth¬ 
ods  of  ten  result  in  a  decrease  in  the  awaienes  of  the  influence 
of  the  specialist’s  decisions  in  his  area  on  other  disciplines.  On 
the  other  hand  it  became  nxrre  and  mofe  tioublesocne  to  ac¬ 
count  strictly  for  all  those  couplings  between  these  subsets 
only  by  parametric  smdies.  In  such  investigations  a  relatively 
small  number  of  principle  parameters  were  varied,  to  find  out 
their  effects  on  the  design  requirements  -  which  were  them¬ 
selves  often  contradictory  -  and  to  improve  the  design. 


Vehicle 

performance 


•  • 


Fig.  1.1:  The  network  of  influences 

With  this  approach,  wmking  in  a  limited  design  space,  the 
engineer  may  achieve  better  results,  but  more  often  it  leads  to  a 
penalty  on  the  design  objectives  to  make  the  initial  concept  fea¬ 
sible. 

A  more  efficient  way  to  integrate  the  diffeier  disciplines 
and  to  balance  their  distribution  in  the  early  design  phases,  is 
the  multidisciplinary  design  optimisation  approach  (MDO). 
Mathematical  optimisation  algorithms  together  with  reliable 
analysis  programmes  and  the  so-called  oprimisation  model 
build  up  a  basis  for  MDO-calculations  with  a  high  rate  of  gener¬ 
ality  and  efficiency.  This  concept  makes  it  possible  to 

•  find  designs  which  meet  all  specified  lequirements  simulta¬ 
neously 

•  achieve  an  optimal  objective  (or  a  combination  of  different 
goals) 

without  time  consuming  manual  and  mote  or  less  intuitive 
search  for  modifications  of  the  initial  design. 

Looking  at  a  typical  data  flow  in  the  structural  design 
phases  of  an  aiictaft  (Fig.  1.2)  the  integrating  effects  of  a  gen¬ 
eral  structural  optimisation  program  can  be  seen.  The  program 
MBB-Lagrange  is  such  a  procedure  which  has  been  developed 
by  MBB  and  sevml  university  institutes  since  1984  [2]. 


In  this  lecture  the  different  parts  and  the  basic  procedure 
of  the  structural  optimisadon  process  are  described. 
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Fig.  1.2:  General  data  flow 


2.  PRACTICAL  ARCHITECTURE  OF 
MULTIDISCIPLINARY  DESIGN 
OPTIMISATION  SOFTWARE 

For  the  treatment  of  opdmisadon  problems  the  "Three- 
Columns-Concept"  [3]  defines  the  pracdcal  architecture  of  an 
opdmisadon  program.  In  the  case  of  structural  design  these 
three  columns  are 

•  Structural  model 

•  Opdmisadon  algorithm 

•  Opdmisadon  model 

The  structural  model  is  the  mathemadcal  descripdon  of 
the  physical  behaviour  of  the  structure,  i.e.  the  necessary 
analysis  procedures  for  calculadng  state  quanddes.  They  are 
often  based  on  finite  element  methods  (FEM)  for  the  struc¬ 
tural  part  and  on  panel  methods  for  the  aerodynamic  calcuia- 
duns  but  other  analysis  methods  can  also  be  applied  (e.g. 
transfer  matrix  procedures  for  special  shell  structures). 

The  opdmisadon  algorithm  is  a  mathemadcal  method  for 
solving  the  general  nonlinear  problem  (NLP). 

minimize  f(x)  (objecdve  function)  (2.1) 

subject  to  g(x)>0  (m^  in-equality  constraints) 

h  (x)  =  0  (m||  equality  constraints) 

Xl^x^Xy  Gower  and  upper  bounds 
for  the  designvaribles  x). 

The  relationship  between  the  structural  model  and  the  op¬ 
dmisadon  algorithm  is  defined  in  the  opdmisation  model, 
which  is  devided  in  the  design  model  and  the  evaluation 
model.  The  design  model  contains  the  transfbnnation  be¬ 
tween  the  mathemadcal  quantities  -  the  design  variables  - 
which  are  processed  by  the  optimisation  algorithm  and  the 
physical  parametera  -  the  structml  variables  -  of  which  the 
optimal  values  have  to  be  determined. 

In  the  evaluation  model,  the  values  for  the  objective  func¬ 


tion  and  to  the  constraints  are  computed  from  the  response 
quanddes  of  the  structural  analysis. 


Fig.  2.1:  Optimization  loop 


Fig.  2.1  shows  the  interaction  of  the  three  columns  in  the 
opdmisadon  process  [4],  First,  the  decision  maker  has  to  de¬ 
scribe  the  structural  and  t^dmisadon  model  for  the  special  de¬ 
sign  problem.  Based  upon  an  initial  design  y^  to  the  structural 
variables,  the  corresponding  itddal  x^^to  the  design  variables 
are  determined.  The  design  model  then  yields  the  variable  sub¬ 
set  of  the  structural  parameters  to  be  opdndsed.  These,  together 
with  the  constant  structural  parameters  (material  constants,  non¬ 
variable  strucnual  parameters),  arc  taken  to  define  a  special  de¬ 
sign  for  which  the  state  variables  are  calculated  by  the  struc¬ 
tural  analysis.  By  means  of  the  evaluation  model  the  objective 
function  and  constraint  values  are  calculated  as  one  pan  of  the 
input  values  to  the  optimisation  algmithm.  In  addition  to  the 
functional  values,  most  optimisation  algorithms  require  the  gra¬ 
dients  of  the  behaviour  functions  with  respect  to  the  design 
variables  which  are  evaluated  by  the  sensitivity  analysis.  If  a 
special  optinnisation  strategy  is  applied,  for  exanqile  a  strategy 
for  solving  a  multicriteria  optimisation  problem,  the  behaviour 
functions  and  their  derivatives  are  transformed  into  correspond¬ 
ing  substitute  values.  Otherwise,  they  are  directly  transferred  to 
the  optimisation  algorithm.  Using  this  information,  the  optimi¬ 
sation  algorithm  calculates  a  new  design  variable  vector  and, 
thereby,  one  obtains  a  closed  optimisation  loop.  If  the  optimal 
design  is  achieved,  which  is  indicated  in  the  optimisation  algo¬ 
rithm  by  breaking-off  criteria,  the  optimisation  process  is  termi¬ 
nated. 

3.  Optimisation  MODEL 

Design  Model 

The  design  trxxlel  describes  the  relationship  between  the 
snnctural  variables  and  the  design  variaUes  determined  by  the 
optimisation  alorithra  If  the  finite  element  method  is  used  to 
structural  analysis  the  following  structural  variaUes  are  possi¬ 
ble. 

•  Sizing  Problems 

Cross  sectional  area  of  elements 
-  Thickness  of  elements 

Laminte  thickness  of  composite  elements 
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Mass  of  concentrated  masses 

•  Geometric  problems 

Fibre  orientation  angles  for  composite  structures 
Coordinates  of  nodes 

Control  parameters  of  parametric  curves 
(e.g.  B-Spline,  Bezier,  polynome, ...) 

Pseudo  loads 

•  Topological  problems 

Arrangement  of  elements  (e.g.  ribs,  spars) 

For  integrated  design  problems,  with  additional  discipline 
analysis  techniques,  such  as  linear  aerodynamics  and  flight 
performance  and  flightmechanics,  new  types  of  variables 
must  be  considered.  Aerodynamic  variables  can  be: 

•  Wing  shape 

Surface  area 
Aspect  ratio 
Taper  ratio 
Sweep  angle 
Profile  shape 

•  Wing  topology 

Arrangement  of  rudders  and  flaps 
Hingelines 

Looking  at  the  flight  performance,  flightmechanics  and 
control  possible  variables  are 

•  Weight 

Gross  weight 
Fuel  weight 
Payload 

•  Wing  load 

•  Mission  parameters 

Range 
Block  times 

•  Thrust  parameters 

•  Fin  volume 

•  Control  parameters 

All  these  different  types  of  physical  variables,  used  in  the 
discipline  analysis  codes,  are  often  not  very  suited  for  a  gen¬ 
eral  mathematical  optimisation  algorithm.  In  order  to  avoid 
numerical  difficulties  and  -  especially  for  the  finite  element 
analysis  -  to  reduce  the  number  of  design  variables,  a  nor¬ 
malization  and  "Linking"  of  variables  is  performed. 

As  an  example,  equation  (3.1)  shows  a  linear  transforma¬ 
tion  between  structural  sizing  and  fiber  orientation  variables 
and  the  corresponding  design  variables: 


with 

t 

vector  of  the  layer  thicknesses  of  all  finite 
elements. 

a 

vector  of  the  layer  angles  of  all  finite  ele¬ 
ments. 

linking  matrices  of  layer  thicknesses  and 
angles. 

x 

design  variable  vector. 

0’  0 

constant  portions  of  die  layer  thicknesses 
and  angles. 

Arbitrary  design  nmdels  can  be  defined  by  the  arrange- 
nnent  of  the  linking  matrices  A^  and  A^^  and  the  vectors  of  con¬ 
stants  and  a^.  It  is  also  possible  to  link  one  comnoon  design 
variable  with  the  structural  variables  of  several  elements  in  or¬ 
der  to  carry  out  a  so-called  "variable  linking".  On  the  other 
hand,  one  structural  variable  depends  at  most  on  one  design 
variable  which  means  that  each  row  of  the  linking  matrices 
contains  at  most  one  coefficient  different  fiom  zero.  The  coef¬ 
ficients  of  the  linking  matrices  and  the  vectors  of  constants  are 
chosen  in  such  a  way  that  the  design  variables  take  on  the  di¬ 
mension  "1”  in  the  design  space  and  are  precisely  ”1”  in  the  in¬ 
itial  design. 

Evaluation  Model 

The  evaluation  model  describes  the  requirements  on  the 
structure  to  be  optimised.  The  special  behaviour,  which  should 
attain  a  minimal  or  maximal  value  in  the  optimal  design,  is 
chosen  as  the  objective  function.  With  aircraft  design  it  is  pri¬ 
marily  important  to  find  a  design  with  a  minimal  structural 
weight.  However,  any  other  state  quantity  (e.g.  costs,  fuel  con¬ 
sumption  etc.)  can  be  considered  as  objective  as  well  if  there 
are  several  objective  functions,  the  problem  has  to  be  solved 
by  multicriteria  optimisation  strategies,  which  are  discussed  in 
15). 

All  nonobjective  requirements  on  the  structural  behaviour 
are  formulated  as  constraints  which  are  normally  upper  and/or 
lower  bounds  on  the  conesponding  state  variables.  For  the  de¬ 
sign  of  aircraft  structures  many  different  types  of  design  re¬ 
quirements  have  to  be  considered  and,  there  still  is  quite  a  lot 
of  work  in  order  to  combine  all  of  the  necessary  analysis  and 
sensitivity  analysis  modules  and  optimisation  nnodules  within 
an  multi^sciplinaiy  t^timisation  system.  The  following  list 
contains  constraint  types  for  the  different  discipline  analyses: 

•  Structural  analysis. 

Thermal  stresses. 

Strength  (failure  safety  factor), 

Di^lacements, 

Stability  fiocal  buckling)  [6, 7], 

Dynamic  quantities  (eigenvalues,  eigenvectors, 
transient-  and  fitequency-response)  [8, 9] 

Manufactoiing  aspects 

•  Including  steady  and  unsteady  aerodynamic  analysis 
methods,  additional  constraints  are  for  example: 

Aeroelastic  efficiencies  [10, 1 1] 

Flutter  speeds  and  danqting 
-  (elastic)  Polar  quantities  (lifi/dtag  ratios) 

(elastic)  Derivatives 
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•  From  flight  perfonnance,  flightmechanics  and  -control 

the  following  constraints  may  arise: 

Manoeuvre  quantities 

>  Roll-  and  turn  rates 

>  Start-  and  landing  perfonnance 

Hinge  moments 

Stability  margins 

Handling  qualities 

The  typical  requirements  for  all  these  types  of  constraints 
is,  that  the  considered  response  quantity  must  be  less  than  an 
allowable  value  (e.g.  stress,  cost)  or  greater  than  a  certain 
limit  (e.g.  flutter  velocity,  roll  rate).  In  the  case  of  the  optimi¬ 
sation  of  a  static  structural  model,  the  number  of  constraints 
can  become  very  high  (e.g.  failure  criteria  in  a  large  multi¬ 
layer  FRP-structure  with  a  lot  of  critical  load  cases  may  lead 
to  100.000  and  mote  constraints).  It  is  clear,  that  dte  treat¬ 
ment  of  such  a  kind  of  problems  is  much  different  to  optimi¬ 
sation  tasks  with  a  few  constraints  only  (e.g.  the  maximum 
turn  rates  of  an  aircraft  or  the  frequencies  of  rigid  or/and  elas¬ 
tic  aircraft  vibrations). 

As  an  example  for  the  mathematical  formulation  of  the 
constraints,  this  latter  mentioned  frequency  requirement  is  de¬ 
scribed  in  the  following: 

gj(x)  :  =  f^  •  f(x)>  0  (3.2) 

or  in  a  normalized  form: 

gj  (X)  :  =  1  -  liii  >  0 

‘max 

This  normalized  representation  has  the  great  advantage  of 
the  independence  of  the  physical  value  of  the  response  quan¬ 
tity  and  guarantees  a  similar  magnitude  for  all  types  of  con¬ 
straints.  Thus  an  improvement  of  the  convergence  behaviour 
of  the  mathematical  optimisation  algorithm  can  be  reached. 

By  means  of  all  these  above  mentioned  constraints  on  the 
state  variables,  many  of  the  most  important  requirements  fw 
the  design  of  aircrafts  can  be  formulated.  All  the  different 
types  of  constraints  and  the  objective  function  -  which  can  be 
defined  by  one  or  more  of  these  constraints  -  form  the  evalu¬ 
ation  model,  and  together  with  the  design  model  they  com¬ 
pletely  describe  the  optimisation  model  and  the  design  task. 

For  practical  applications  it  can  not  be  expected  and  even 
it  is  not  desired  that  there  will  be  one  computer  program  only 
for  the  optimal  design  of  aircraft.  A  much  better  solution  for 
the  multidisciplinary  design  optimisation  of  a  complex,  inter¬ 
nally  coupled  system  behaviour  is  the  separate  evaluation  of 
the  individual  discipline  analysis  and  the  partial  sensitivity 
analysis  with  a  well  organized  exchange  of  input  and  ouqrut 
data.  But  most  important  is  an  efficient  method  for  calculat¬ 
ing  the  coupled  system  sensitivities. 

Such  a  formulation  is  presented  in  [13]  and  shortly  de¬ 
scribed  in  the  second  part  of  diis  lecture  by  iq^lying  it  to  the 
integrated  design  of  a  fin. 


4.  THE  OPTIMISATION  ALGORITHM 

Another  column  in  the  Three-Column-Concept”  repre¬ 
sents  the  cqrtiniisation  algarithms.  In  the  previous  chatter  the 
different  types  of  problems  in  the  multidi«aplinary  design  op¬ 
timisation  pnxess  were  shown.  Many  practical  applications  in 
the  last  dMade  have  proved,  that  it  is  necessary  to  inovide 
several  different  tqmmisation  strategies  and  algorithm  to  get 
reliable  solutions,  because  there  is  no  known  single  method 
which  is  ad^ted  to  every  type  of  proUenL 

To  understand  the  solution  process  for  the  NLP-problem 
fonnulated  in  equaticm  (2.1),  it  is  necessary  to  formulate  the 
required  optimality  ctmditions  (Kuhn-Tucker-Conditions); 

V,  L  (X*,  X*)  =  0  <=>  7,f(x*)  =  ^(x*) 


gj  (X*)  =  0 

X*j  >  0 

(4.1) 

where 

II 

is  the  Lagrangian  function 

X* 

the  c^mmal  solution  vector 

the  Langrangian  multiplier  in 
the  optimum 

and 

^x 

the  gradient  with  respect  to  x. 

That  means,  that  in  the  optimum  the  gradient  of  the  objec¬ 
tive  function  is  a  non-negative  linear  combination  of  the  gra¬ 
dient  of  the  so-called  "active  constraints".  The  determination 
of  these  active  constraints  is  one  of  the  main  problems  of  all 
(^timisation  procedures.  The  "less  active"  constraints  have 
less  influence  on  the  current  detign  change  and  are  therefore 
temporarily  neglected.  Suitable  deletion  of  these  constraints 
accelerates  the  optimisation  process,  but  it  is  not  easy  to  man¬ 
age. 

To  find  the  optimal  solution  vector  x*,  most  of  the  mathe¬ 
matical  proramming  algorithms  uses  the  following  iterative 
formulation: 


where  is  the  downhill  search  direction  and  the  step 
size,  is  a  positive  scalar,  which  minimizes  a  function  F  in 
the  direction  of  sousing  a  one-dimensional  line  search,  that 
means: 

F(x'^')  =  min  [F(x''-t- a''s'')]  (4.3) 

a 

The  formulation  of  F  depends  on  the  optimisation  method 
and  is  explained  somewhat  later.  The  calculation  of  the  step 
size  o'  is  a  relatively  simple  matter,  which  however  requires 
the  evaluation  of  the  structural  model  and  must  therefore  be 
carried  out  very  effeciently. 

Without  going  too  much  into  detail,  a  classification  of 
mathematical  programming  methods  is  given  below  [14]: 

•  Transformation  methods 
Penally  functions 
Barrier  functions 
Method  of  multipliers 
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For  these  strategies,  the  function  F  is  the  original  objec¬ 
tive  function  f  (x),  augmented  by  a  weighted  penalty  term, 
which  summarizes  the  constraint  values  g  (x). 

For  the  inverse  barrier  method  (IBF)  for  example,  the 
transformed  problem  can  be  written  as 

F  (x'')  =  f(x'')-i-r''j:  gy‘(x'')  (4.4) 

where  the  penalty  parameter  r'*  is  updated  after  each  iteration 
step  by 

j.v+1  3-  jVp^  0<c<l  (4.5) 

These  strategies  are  very  reliable,  but  they  need  a  lot  of 
function  and  gradient  evaluations  and  arc  very  "expensive". 

•  Primal  methods 

Indirect  methods 

Sequential  linear  programming  (SLP) 

Sequential  quadratic  programming  (SQP) 

These  methods  solve  a  sequence  of  linearized  or  quad¬ 
ratic  subproblems.  In  the  sequential  linear  programming 
methods  the  nonlinear  functions  f(x)  and  g(x)  are  expanded  in 
Taylor  series  considering  only  the  linear  terms  (Fig.  4. 1 ).  This 
linearized  problem  can  then  be  efficiently  solved  by  using  the 
simplex  algorithnL  For  problems  with  many  active  constraints 
this  method  works  very  efficiendy.  (Normally  less  than  10  it¬ 
erations  ate  needed!).  But  if  there  are  only  a  few  active  con¬ 
straints,  convergence  can  worsen.  In  the  case  of  highly  non¬ 
linear  problems  (e.g.  buckling,  structural  dynamics),  the 
method  tends  to  fail  because  of  the  tough  approximation  (= 
linearization)  of  the  original  problem. 


Fig.  4.1:  Sequential  linearization  SLP 


For  the  SQP-methods,  the  quadratic  subproUems  result 
from  a  second  order  approximation  of  the  Lagiangian  func¬ 
tion  L(xA)  and  a  linearization  of  the  constraint  functions  g(x). 
The  search  direction  a  of  equation  (4.2)  is  then  found  by  solv¬ 
ing  the  following  quadratic  subproUem  [15]; 

ntin.  [  !^4  8^  f  (x'')  1  (4.6) 

7,gj(x'')'^a  gj(*'')2  0 

where  is  an  qrproximation  of  the  Hessian  matrix  of  the 
Lagranghui  function  for  Ae  v-th  iteration  step. 


These  methods  ate  very  general  and  robust.  They  can  be 
used  for  a  wide  range  of  problems,  independent  of  the  ratio  of 
active  constraints  and  design  variables.  The  accuracy  of  the 
optimal  result  is  very  good.  If  the  starting  point  is  far  from  the 
c^timal  solution,  it  might  happen,  that  a  lot  of  iteration  steps 
(=  function  and  gradient  evaluations)  are  needed  (ttxrte  than 
20).  For  this  reason,  the  best  efficiency  is  shown  for  medium 
size  problems. 

•  Direct  methods 

-  Gradient  projection  method  (GPM) 

-  Generalized  reduced  gradients  (GRG) 

Method  of  feasible  directions  (MFD) 


One  principle  idea  of  the  GRG-method  is  the  transforma¬ 
tion  of  the  m^  inequality  constraints  g(x)  of  the  original  prob¬ 
lem  (2.1)  into  equality  constraints  by  introducing  additional 
variables.  By  this  means  the  c^timisation  process  is  working 
in  the  feasible  domain.  This  leads  to  the  following  ttKxiiried 
optimisation  problem: 

min.  f  (x) 

s.  t.  h  (x)  =  0  (4.7) 

;  i=l...n 

0<  X-  <  oo  ;  i  =  n-t-l,n-fmg 

This  system  of  (m^  +  m^)  equations  with  (n  ■(■  m^)  un¬ 
knowns  gives  a  solution  for  +  m„  so-called  basis  vanables 
which  depends  on  (n  -  mj,)  nonbasis  variables.  A  clever  sepa¬ 
ration  technique  for  these  two  types  of  variables  and  a  lineari¬ 
zation  of  the  constraint  functions  b(x)  results  in  a  linearized 
objective  function  f|^(x),  with  a  reduced  set  of  variables,  which 
depend  only  from  the  (n  -  m,,)  non-basis  variables  x.  The 
search  direction  s''  for  this  smaller  problem  can  be  found  for 
example  by  using  the  negative  reduced  gradient  of  the  objec¬ 
tive  function 


(4.8) 


or  a  modified  direction,  which  take  into  account  second  order 
infonnations. 

An  important  improvement  of  the  efficiency  of  the  GRG- 
method  is  the  use  of  a  SQP-seatch  direction  that  means  the  so¬ 
lution  of  a  quadratic  subproblem  [16],  This  hybrid  SQP-GRG- 
Algorithm  can  reduce  considerably  the  number  of  furtetion 
calcs  and  gradient  evaluations. 

•  Dual  concepts 

The  priiKiples  of  the  dual  formulation  ate  summarized  in 
the  following  [17];  The  solution  of  the  primal  problem  (equa¬ 
tion  2.1)  can  be  obtained  by  a  "Min-max"  two  phase  proce¬ 
dure: 

maximize  1  (A.) 

subject  to  Xj  2  0 

where  the  dual  function  1(X)  which  depends  only  on  the 
Langrangian  multipliers,  result  from  minimizing  the  Lagran- 
gian  L  (x,  X)  over  the  allowaUe  primal  variables; 

I(X)  »  min  L(x,X)  (4.9) 

Very  important  for  the  practical  application  of  this  method 
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is  the  approximation  concept  of  the  objective  and  the  con¬ 
straint  functions.  A  so-called  "convex  linearization"  (18)  for 
example,  leads  to  a  sequence  of  convex  and  separable  sub- 
problems  with  a  likewise  separable  Lagrangian  function, 
which  can  be  solved  easily  by  one-dimensional  minimizing 
methods; 

min  L:  (x.)  =  3;  x,  (4.10) 

S.L  Xlj<Xi<X„i  ’'i 

where  the  coefficients 


f; 


and  bj  =  ^djjXj 


depend  only  upon  the  dual  variables  Xj.  The  coefficients 
fj,  djj  and  C;j  results  from  a  mixed  approximation  of  the  objec¬ 
tive  and  the  constraint  functions,  where  the  fj  represents  the 
first  derivatives  of  the  objective  function  and  the  d-  denote 
the  first  derivatives  of  the  constraint  functions  with  respect  to 
the  design  variables  x,  (i.e  the  components  of  the  gradients). 
The  Cjj  are  the  first  derivatives  of  the  constraint  functions 
with  respect  to  the  reciprocal  variable  z,  =  l/x^.  (The  type  of 
the  constraint  approximation  (direct  or  reciprocal)  can  be  de¬ 
cided  by  the  sign  of  the  derivatives  for  example). 

Practical  applications  of  the  (X)NLIN-algorithm  for 
structural  optimisation  problems  have  shown  a  very  efficient 
convergence  behaviour  for  sizing  as  well  as  for  shape  design 
tasks. 

Besides  these  mathemtical  programming  methods  there 
exists  another  approach  to  solve  the  structural  optimisation 
problem  -  the  optimality  criteria  procedure  (19). 

For  this  formulation  the  stationary  conditions  of  the  La¬ 
grangian  function  in  the  optimum  (Kuhn-Tucker-conditions, 
eqs.  4.1). 


=  ZX:  .51 
5x  J  J  Ox 

is  written  in  the  form 


?eyX,=  1 


where  eji  is  the  ratio  of  the  first  derivatives  of  the  con¬ 
straints  and  me  objective  function.  This  set  of  equations  can 
be  solved  easily  for  the  unknown  Lagrangian  multiplies  Xj 
with  a  kind  of  separability  assumption  for  the  active  con¬ 
straints,  which  leads  to  an  estimate  of  the  Lagrangian  multi¬ 
pliers: 


where  a  is  relaxation  parameter  and  can  be  seen  as  a  step 
size  parameter  (e.g.  a  =  2  assures  a  reasonable  rate  of  con¬ 
vergence).  The  optimality  criteria  approach  lead  with  a  reltive 
small  amount  of  computing  effon  to  a  solution,  almost  re¬ 
gardless  of  me  number  of  variables.  This  is  in  general  how¬ 
ever,  not  an  optimal  design,  especially  if  mere  exists  not  only 
one  dominant  type  of  constraint  in  the  optimum,  but  a  variety 
of  different  constraint  types,  which  is  often  me  case  in  mutli- 
disciplinary  design  optimization  problems.  Therefore  it  is 
necessary  to  have  a  very  good  understanding  of  me  physics  of 
the  problem  to  decide,  if  an  optimality  criteria  inemod  can  be 
used. 


5.  DISCIPLINE  ANALYSIS  AND  SENSI¬ 
TIVITY  ANALYSIS 


The  task  of  discipline  analysis  is  to  calculate  me  state 
quantities  of  the  strucnire  required  to  determine  me  constraint 
and  objective  values  defined  in  me  evaluation  model.  As 
mentioned  before,  most  optimization  algorimms  do  not  only 
require  me  functional  values  of  me  behaviour  functions  but 
also  meir  sensitivities  wim  respect  to  me  design  variables. 
The  calculation  of  mese  sensitivities  can  be  carried  out  ana¬ 
lytically  as  well  as  numerically  by  means  of  simple  differen¬ 
tial  quotients.  Since  in  me  aircraft  design  one  has  often  large 
scale  design  problems  wim  sometimes  several  hundreds  of 
design  variables,  it  is  necessary,  fen-  me  sake  of  calculation  ef¬ 
fort  and  economy,  to  determine  me  sensitivities  by  me  ana¬ 
lytical  differentation  of  all  descriptive  equations  as  far  as  pos¬ 
sible. 

In  me  following  a  brief  survey  of  different  type  of  analy- 
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s  and  sensitivity  calculations  is  given. 


The  structural  analysis  is  based  on  me  finite  element 
rnemod.  This  is  a  well-known  reliable  and  very  general  way 
of  modelling  bom  me  static  and  me  dynamic  behaviour  of 
structures.  It  is  possible  to  treat  h<»ix>geneous  materials  wim 
isomropic,  ormotropic  and  anisotropic  pre^rties  as  well  as 
composite  materials.  (For  special  types  of  structures  it  can 
make  sense  to  use  omer  -  often  very  efficient  -  methods  to  de¬ 
scribe  me  re^;>onse  quantities,  e.g.  Kitchhoff  plate  theoiy  for 
min  wing  structures  [20]  or  transfer  matrix  procedures  for  cy- 
lindric  shell  structures  (21). 


StotkoroMenis 

For  linear  elastic  structures  wim  static  loads,  the  funda¬ 
mental  stiffness  equation  describes  the  structural  response: 


Now  an  iterative  resizing  algorimm  can  be  derived  by 
multiplying  bom  sides  of  equation  (4.11)  by  x°  and  taking 
me  a-m  root: 


K(x)  u(x)  =  p  (tt,x) 


K  (x)  stiffness  matrix 


u  (x)  displacement  vector 
p  (u,x)  load  vector 


is  composed  of  a  pan  that  depends  on  the  solution  and  one  that 
is  independent  of  it.  Thus  equation  (5.1)  can  be  written  as: 

K(x)u(x)  =  p„+p(u)  (5.4) 


With  the  displacement  vector  u  (x)  and  the  design  vari¬ 
ables  X,  the  stress  and  the  strains  in  the  structure  can  be  calcu¬ 
lated.  as 

Besides  these  strength  quantities,  for  lightweight  aero¬ 
space  structures  exists  the  important  problem  of  local  and 
global  instability,  than  means  large  deformations  of  the  struc¬ 
ture.  Two  main  concepts  for  solving  the  buckling  mechanism 
are  shortly  described  in  the  following  (Fig.  5.1). 


Assuming  linear  aerodynamics  the  load  can  be  expressed 
Pq  +  p  (u)  =  qTjL  F  •  A  (w^j  -I-  w^)  -t  (5.5) 

with 

q  =  -£  v^  Dynamic  pressure 


Fig.  5.1:  Concept  for  stability  calculation 


p  Air  density 

V  Air  speed 

Tjl  Transfoimation  of  aerodynamic  panel 

forces  into  finite  element  mode  forces 

F  Aerodynamic  panel  surfaces 

A  Aerodynamic  influence  matrix 

w^  Angle  of  attack  of  the  panels  for  the 

rigid  structure 

Wp  Angle  of  attack  of  the  panels  due  to  elastic 

deformations  of  the  structure 

Mass  loads 


One  possibility  is  the  formulation  of  an  eigenvalue  pro- 
belm: 

(K  (X)  +  X  •  Kg  (x)|  =  0  (5.2) 

where 

Kg  (X)  is  the  geometrical  stiffness  matrix  the 

lowest  eigenvalue,  which  defines  the  cri¬ 
tical  load  by  Pj,|.  =  X  •  p 

and 

u  (x)  the  eigenvector,  which  represents  the 

buckling  shape. 

For  local  instability  problems  (buckling  of  bars  and 
shells)  the  critical  loads  and  stresses  can  often  be  calculated 
with  special  stability  equations  (e.g.  the  well-known  Euler 
equations  for  bars)  (6,7). 

In  the  case  of  a  two-dimensional  loading,  considering  ten¬ 
sion/compression  and  shear  forces  with  the  following  formula 
the  margin  of  safety  can  be  defined: 

|g  (x)|  ^  T(x)^  ^  j  (5.3) 

Ocr  ^cr 

The  calculation  of  the  critical  stresses  o^,  can  be  done  by 
solving  analytically  the  differential  equations  of  the  plate- 
theory  with  special  material  and  geometrical  assumptions. 

In  the  case  of  aeroelastic  problems,  the  load  vector  p(ii,x) 
depends  on  the  deformation  of  the  structure.  The  load  vector 


Expressing  the  elastic  pan  by  the  finite  element  node  dis¬ 
placements,  the  aerodynamic  forces  acting  at  the  nodal  points 
of  the  FE-mesh  become 

Pa  =  q  Tsl  F  A  Wq  q  Tjl  F  A  Tlc  Tjc^u  (x)  (5.6) 

where  the  transformation 

T^f-  relates  the  panel  comer  displacements  with  the 
equivalent  panel  angle  of  attack 

and 

Tjp  the  panel  comer  displacements  with  the  finite 
element  node  displacement 

With  the  abbrevation 

the  aeroelastic  equilibrium  equation  can  be  written  in  the  gen¬ 
era)  form 


(K-C)u(x)  =  Pq  (5.7) 

If  the  difference  (K  -  C)  is  regular,  that  means  non¬ 
singular,  this  equation  can  be  solved.  Because  of  the  non¬ 
symmetry  of  the  aerodynamic  influence  matrix,  which  is  part 
of  the  matrix  C,  the  numerical  effon  fw  a  direct  solution  of  the 
equation  (5.7)  became  very  high,  already  for  minor  problems. 

For  that  reason  an  iterative  solution  procedure  was  devel¬ 
oped,  where  an  additional  relaxation  process  is  introduced  to 


improve  convergence  110]: 

K  =  (0 C  +  (l-co)  K  u®  +  0)  Pp  (5.8) 

The  convergence  of  the  iteration  strongly  depends  on  the 
dominant  eigenvalues  X,  of  the  corr-sponding  eigenvalue 
problem: 

(toC  +  (l-o))K-X,(a))K)w  =  n  ,  (5.9) 

with  the  eigenvalue  transformation 


lion  with  the  tight  side  of  equation  (5.8)  using  the  vector  u*'* 
of  the  proceeding  iteration  step.  Practical  appl’cations  have 
shown  a  good  convergence  behaviour  of  this  solution  process. 
With  the  optimal  value  of  the  relaxation  parameter  to,  the  dis¬ 
placement  vector  u(x),  due  to  aerodynamic  and  mass  forces 
can  be  determined.  The  solution  can  then  be  used  fw  the  com¬ 
putation  of  a  so-called  static  aeroelastic  efficiency  of  a  struc¬ 
ture.  These  factors  describe  the  influence  of  the  elastic  struc¬ 
ture  on  the  aerodynamic  forces  and  moments  an  are  usually 
expressed  in  the  form 

total  load 

^  rigid  load  (5.15) 


Xj(a))  =  X(o  -  CO  -H  1  (5.10)  This  ratio  is  normally  less  than  one  and  depends  on  the 

dynamic  pressure  q  acting  on  the  structure.  From  equations 
where  X  is  either  the  maximal  or  the  minimal  eigenvalue  of  ^  aeroelastic  efficiency  can  be  obtained 

the  original  problem  (eq.  5.7) 


(C  -  X  K)  w 


(5.11) 


_  -  [pn  -I-  C  u(x)] 

s^Po 


(5.16) 


which  can  be  computed  by  a  simple  v.  Mises-lteration: 

Ax®''  =  K-'CAx®  (5.12) 


An  approach  for  the  optimal  relaxation  parameter  co  can 
than  be  found  by  the  mean  of  the  minimal  and  maximal  ei¬ 
genvalues  from  the  following  simple  geometric  relation, 
which  is  also  indicated  in  Fig.  5.2: 


co„„  = 


opi  2  -  X  -  X 

^min 


(5.13) 


In  most  cases  K  represents  the  symmetric  and  banded  stiff¬ 
ness  matrix  of  the  finite  element  model  of  the  structure  (e.g.  a 
wing)  which  can  be  decomposed  by  a  Cholesky  facturization 
as  follows: 


K  =  L 


(5.14) 


where  L  is  a  lower  triangular  matrix.  This  fact  is  very  impor¬ 
tant  for  the  efficiency  of  the  methcxl,  because  one  iteraion 
step,  that  means  a  bener  approximation  of  the  solution  vector 
u'"'*'  %  requires  only  one  forward  and  one  backward  substitu- 


with  s  as  a  vector  for  summing  up  the  forces  or  moments  of 
those  aerodynamic  panels  which  contributes  to  the  efficiency 
(e.g.  the  panels  on  a  control  surface). 

By  using  only  unit  cases  for  the  angle  of  attack  of  the 
panels  i.g.  =  1,  it  is  possible  to  compute  all  a-dependent 
derivatives  of  the  elastic  complete  aircraft  (e.  6-  ^-La’  ^ma 
and  by  the  same  way  derivatives  with  respect  to  the  side¬ 
slip  angle  P  (e.g.  •.-Yp:  c®n  be  calculated,  too. 

These  method  is  also  applicable  to  the  determination  of 
the  derivatives  which  depend  on  the  rotational  degrees  of 
freedom  of  the  airplane,  i.e.  roll  velocity  p,  pitch  velocity  q 
and  yaw  velocity  r.  For  these  cases  the  distribution  of  angle  of 
attack  a  and  sideslip  p  depend  on  the  distance  of  the  panels  to 
the  corresponding  axis,  respectively  to  the  center  of  gravity 
and  on  the  flight  velocing  v.  With  that,  the  angle  of  attack  Oj 
of  a  panel  i  due  to  pitching  for  example,  can  be  written  as: 

aj  =  q’‘i/v  (5,16) 

with  Xj  as  distance  fiom  the  pitch  axis. 

To  calculate  these  clastic  derivatives  of  the  complete  air¬ 
craft,  the  decomposition  of  the  stiffness  matrix  has  to  be  done 
(eq.  5.14),  which  is  only  possible  if  the  matrix  K  is  positiv  de- 
finit  That  means,  that  the  airplane  either  has  to  be  supported 
or  the  rigid  body  degrees  of  freedom  must  be  eliminated  by 
special  expensive  and  time  consuming  transformations  [22],  It 
can  be  shown  [23],  that  a  statical  determined  suppmt  of  the 
aircraft  (e.g.  close  by  the  center  of  gravity)  gives  correct  re¬ 
sults  for  total  aircraft  loads  using  the  unit  case  method,  if  trim 
conditions  are  considered.  That  means,  that  the  total  sum  of 
forces  and  moments  due  to  aerodynamics  and  masses,  which 
act  on  the  airplane,  has  to  be  zero.  These  trim  conditions  can 
be  written  in  the  following  short  form  [12],  [23]. 

z'''r  =  z’''H'I'-i-ngZ'''h„  (5.17) 


with 
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r 


H 


and 


total  sum  of  external  loads 
unit  aerodynamic  loads 
mass  loads 
load  factor 

marix,  which  contains  the  trim  conditions 
(=  sum  of  forces  and  moments) 


and 


(5.23) 


^  =  -  (K  -  C)''  ^  u 
5xj  6Xj 

and  using  an  auxiliary  vector  d,  which  can  be  calculated  by  the 
same  iterative  solution  procedure  as  shown  in  equation  (5.8) 


V 


the  vector  of  the  unknown  factors  for  the 
trim  parameters  (e.g.  rudder  deflection,  an¬ 
gle  of  attack  for  a  steady  two  degree  of 
freedom  longitudinal  case). 


h^(K-C)'  =  d 

the  derivative  can  be  calculated  finally  by 


The  total  elastic  deformation  u(x)  is  finally  achieved  by  5g  _  __1_  6K  ^ 

multiplying  the  deformations  iiuj.(x)  resulting  from  the  unit  "Sxj  6xj 

cases  with  the  scaling  factors 


(5.24) 


(5.25) 


u(x)  =  u^(x)  'Ffx)  -t-  ng  u„(x)  (5.18) 

where  u^(x)  is  the  deformation  due  to  mass  load  for  load  fac¬ 
tor  one. 

Using  the  displacements  determined  by  the  global  static 
structural  analysis,  the  strains  and  stresses  can  now  be  calcu¬ 
lated.  Especially  for  fiber  compostite  structures  the  safety 
against  material  failure  is  usually  checked  by  means  of  vari¬ 
ous  failure  criterias,  e.g.  according  to  Tsai-Wu,  Tsi-Hill, 
Hoffmann  and  others  [24], 

All  these  in  the  foregoing  sections  explained  state  vari¬ 
ables,  which  will  be  generally  denoted  by  the  vector  r  in  the 
following,  depend  on  design  variables  with  an  explicit  de¬ 
pendency  on  the  equation  parameters  and  on  implicit  depend¬ 
ency  on  the  structural  defoimation  n.  Therefore  the  corre¬ 
sponding  constraints  are  formulated  as 


The  derivtive  of  the  stiffness  marix  K  can  be  calculated 
analytically  for  thickness  and  fiber  orientation  as  design  vari¬ 
ables  by  using  the  finite  clement  formulation  [25].  For  general 
geometry  variables  the  sensitivity  can  he  achieved  numeri¬ 
cally; 

6K  _  K  (x  +  EXiCj)  -  K  (x) 

flXj  £  Xj  (5.26) 

with 

Cj  unit  vector 

c  small  real  number 


The  computing  effort  for  the  additional  stiffness  matrices 
can  be  essentially  reduced,  if  only  the  terms  effected  by  Xj  are 
calculated  anew. 


g  =  g[r(x,u)]  (5.19) 

and  the  derivatives  of  the  constraint  vector  g  with  respect  to 
the  design  variables  as  it  is  needed  for  the  optimisation  algo¬ 
rithm  can  then  be  achieved  by  using  the  chain  rule: 

"Sx  6r  \8x  Ou  0x1  (5.20) 

The  derivative  of  the  constraint  vector  g  with  respect  to 
the  state  variables  r  depends  only  on  the  applied  discipline 
analysis  and  is  determined  by  the  evaluation  model.  In  the 
following  the  solution  method  is  shortly  explained  for  the 
state  variable  aeroelastic  efficiency  t|,  where 

g=2L  -  1  >0  (5.21) 

n 

'nun 

Differentation  with  respect  to  design  variable  x.  leads  to 

^  ^  _1_  ^  ^  (5.22) 

I  'mm  I 

With  equations  (5.7,5.15) 

-511= 

Ou 


Dynamic  Problems 

The  general  equation  of  motion  describing  time-dependent 
system  deformations  reads  as  follows: 

Mii-t-Du-fKu  =  p(t)  (5.27) 

with 

M  Mass  matrix 
D  Damping  matrix 

t  Time 

In  the  case  of  undamped  eigenvibrations,  equation  (5.27) 
can  be  transformed  into  a  real  eigenvalue  problem  using  a  har¬ 
monic  approach  for  the  displacements 

(K-(OjM)yj=0  (5.28) 

where 

(Uj  is  the  j-th  natural  frequency 
and  yj  the  j-th  eigenvector 

The  constraints  on  natural  frequencies  usually  consist  in 
imposing  lower  or  upper  limits.  With  the  normalized  formula¬ 
tion  (eq.  3.2,  5.19)  the  sensitivity  of  problems  with  frequency 
constraint  is  given  by 
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Sgj  =  Ogj  1 
5)i~  2to. 


CK 

8x 

1 


r^\ 

axilj 


(5.29) 


The  derivative  of  the  mass  matrix  M  can  be  calculated  in 
the  same  way  as  described  earlier  for  the  stiffness  matrix  K. 

A  more  advanced  type  of  constraints  from  dynamics  is 
given  by  frequency  or  transient  response  problems.  These 
may  occur  when  a  structure  is  loaded  harmonically  or  by  a 
time  dependent  load  and  it  is  required  that  the  displacement, 
velocity  or  acceleration  at  certain  points  of  the  structure  must 
not  exceed  prescribed  values. 

To  get  the  structural  response  quantities  in  the  time- 
domaine  the  second  order  differential  equation  (5.27)  has  to 
be  solved.  This  can  be  done  for  example  by  the  method  of 
Newmark. 

The  displacement  vector  of  a  harmonically  loaded  struc¬ 
ture  (excitation  frequency  fl)  is  computed  by 

(-Q^M  -t-  iQD  -K  K)  u  =  p(fl)  (5.30) 


The  solution  of  these  equations  (5.27,  5.30)  is  quite  ex¬ 
pensive.  To  reduce  the  computational  effort,  normally  a  trans¬ 
formation  is  introduced: 

u  =  Tq  (5.31) 

The  transformation  matrix  T  is  a  n  x  m  -  Matrix  (m 
«n)  which  reduces  the  original  system  drastically.  If  T  con¬ 
tains  eigenvectors  of  the  structure  (normally  the  lowest  one), 
equation  (5.31)  represents  the  transformation  to  modal  coor¬ 
dinates.  In  recent  years  the  transformation  to  Lanczos  coordi¬ 
nates  have  been  proposed  in  structural  dynamics.  This 
method  was  shown  to  be  especially  promising  because  it  is  a 
load-depenent  transformtion  which  approximates  the  influ¬ 
ence  of  higher  modes  as  well. 


Modaltraasform»tioo 
5  EigMvektoreo 


Mod*Itr&iufonn4iioa 
20  Eisenvektorea 


L&nczoC'TransformAtiofi 
5  LiAflczocvektoren 


LAQcz0t*Tr4asform4tioo 
20  Luczosvektoreo 


Fig.  5.3  shows  the  results  for  the  acceleration  of  a  point  of 
a  cantilever  plate  in  the  time  domaine.  In  this  example  five 
and  twenty  transformation  vectors  for  the  modal  respectively 
the  Lanczos  transformation  are  used 

It  can  be  seen  that  only  five  vectors  on  the  Lanczos  case 
are  necessary  to  get  good  results  (gj.  Formulation  of  con¬ 
straints  on  transient  and  frequency  response  quantities  and  the 
corresponding  sensitivity  analysis  is  quite  complicated  and 
out  of  scope  of  this  lecture  and  well  described  in  |8,9J. 

A  special  kind  of  barmonically  loading,  however,  must  be 
mentioned  in  this  context  The  phenomena  of  self-exciting  vi¬ 
brations  of  elastic  structures  in  a  flow  field,  which  is  a  dy¬ 
namic  stability  problem  and  designated  as  flutter.  Due  to  the 
interaction  of  the  aerodynamic  forces,  the  elastic  forces  and 
the  inertia  forces  with  the  structural  deformation,  there  is  an 
exchange  of  kinetic  energy  of  the  air  flow  with  the  elastic  and 
kinetic  energy  of  the  structure 

At  the  boundary  between  damping  and  excitation  there  is 
no  energy  exchange,  which  means  that  small  disturbances 
lead  to  harmonic  vibrations.  Depending  on  the  stiffness  and 
mass  distribution  of  a  structure,  such  a  critical  case  occurs 
when  certain  combinations  of  flow  velocity  and  Mach  number 
are  given.  The  corresponding  critical  flow  velocity  is  called 
flutter  speed.  Since  no  flutter  case  can  be  admined  in  the 
whole  mission  ra.nge  of  an  aircraft  it  must  be  required  that  the 
smallest  fluner  speed  does  not  fall  shon  of  a  certain  limit 
given  by  the  maximal  flight  speed  plus  safety  increase  (15% 
safety  increase  for  military  aircrafts,  20%  for  civial  aircrafts). 
The  maximal  flight  speed  can  be  taken  from  the  so-called 
flight  envelope,  which  depicts  the  mission  range  of  an  air¬ 
craft. 

According  to  (2. 1 )  the  fluner  constraint  can  be  formulated 
as 

g  =  ''F/v„«„-l  (5.32) 

with  Vp  flutter  speed 

^max  tnaxinium  flight  speed 

For  the  determination  of  the  flutter  point,  that  means  the 
calculation  of  the  critical  flow  velocity  Vp,  harmonic  aerody¬ 
namic  forces  which  depend  on  the  harmonic  deflection  u  of 
the  structure  are  introduced: 

p(t,0))  =  C(o>,Ma)  u  e'“*  (5.33) 


where  C  contains  the  conqrlex  aerodynamic  influence 
matrix  and  the  transformation  of  loads  from  the  aerodynamic 
into  the  finite  element  mesh  similar  as  described  far  the  static 
aeroelastic  (eq.  5.6)  and  the  dynannic  pressure. 

The  aerodynamic  influence  matrix  is  a  fully  occupied, 
non-symmetric  complex  matrix  depending  on  the  Mach  num¬ 
ber  and  the  reduced  frequency  k.  For  constant  altitudes  the 
dependency  on  the  Mach  Number  can  be  transformed  into  a 
dependency  on  the  airspeed.  With  (5.27, 5.33)  and  by  neglect¬ 
ing  material  damping  the  flutter  analysis  equation  can  be  writ¬ 
ten  as: 

[Jf-C(vJt)  +  xW)A#]q(vJc)  (5.34) 


Fig.  S3i  Comparison  of  the  acceleration  of 
a  cantilever  plate 


This  equation  contains  a  system  reduction  according  to 


2-11 


(5.31)  to  reduce  the  numerical  effort, with 


K  =  T^KT 

Generalized  stiffness  matrix 

M  =  T^MT 

Generalized  mass  matrix 

C  =  T^CT 

Generalized  aerodynamic  load  matrix 

1 

q 

generalized  coordinates,  right  hand 

eigenvector 

The  complex  eigenvalue  problem  of  equation  (5.34)  can 
be  solved  for  example  with  a  QR-algorithm  or  -  if  an  initial 
solution  is  known  -  with  the  very  efficient  perturbation 
method  by  Wittmeyer  [26]. 

The  resulting  complex  eigenvalues  of  the  flutter  equation 
depend  on  the  reduced  frequency  k  and  the  velocity  v: 


which  the  real  part  of  the  eigenvalue  vanishes,  that  means 
an  undamped,  harmonic  motion  takes  places  (Fig.  5.4)  with 
the  bequency  Xp”. 

Using  this  flutter  velocity  in  the  constraint  equation  (5.32), 
the  derivative  of  the  flutter  constraint  with  respect  to  design 
variables,  can  be  achieved. 

At  first  the  differentation  of  the  flutter  equation  (5.34) 
leads  to  the  sensitivity  of  the  eigenvalue  in  the  flutter  point: 


^  =  -_L  pVM  -  .5£  -h 

6Xj  2Xp  y6Xj  6x.  fiXj/ 


(5.38) 


with 

p  left  hand  eigenvector 


and  the  normalization 


Xp(v,k)  =  Xp’(vjc) -t- i  Xp”(v,k)  (5.35) 


where 

is  the  real  part  or  the  damping 

The  sensitivity  equation  of  the  flutter  speed  itself  is  ob¬ 
tained  by  re-arrangement  and  differentation  of  the  definition 
equation  of  the  reduced  frequency: 

the  imaginary  part  or  the  eigenfrequency. 

With  the  notation  of  (5.35)  the  definition  equation  of  the 
reduced  frequency  k  of  the  aerodynamic  matrix,  is  given  by 


6vp  _  1  6Xp  1  Xp  6k 
6Xj  k  6Xj  k^  6Xj 


(5.39) 


k  =  (5.36) 

V 

or 

XA”(vJc)  =  -ili- 

That  means,  that  valid  points  for  a  flutter  curve  are  only 
those,  where  the  imaginary  part  of  a  solution  Xp,  i.g.  the  fre¬ 
quency  Xp”,  corresponds  to  the  frequency  X^”  of  the  oscil¬ 
lating  airload.  This  requirement  can  be  formulated  by  the  fol¬ 
lowing  intersection  condition 

X”  =  X^”  =  Xp”  =  X”[v(k)j  (5.37) 


Fig.  5.4:  Evaluation  of  flutter  curves 

A  flutter  point,  bnally.  is  found  for  a  velocity  Vp  in 


With  these  basic  equations  (5.38,  5.39),  and  the  additional 
condition,  that  in  the  flutter  point,  the  real  pan  of  the  eigen¬ 
value,  i.g.  the  damping  vanishes  and  finally  taking  into  ac¬ 
count  the  dependency  of  the  aerodynamic  influence  matrix 
from  the  reduced  frequency  k  and  the  Mach  number  Ma,  the 
derivation  of  the  flutter  constraint  can  be  obtained. 

It  should  be  mentioned,  that  the  differentiation  of  the 
transformed  matrices  K,  M  and  C  with  respect  to  the  design 
variables  includes  terms,  which  depend  on  the  derivative  of 
the  transformation  matrix  T.  In  the  case  of  a  iixxlal  transfor¬ 
mation,  T  contains  a  number  of  eigenvectors  of  the  undamped, 
homogeneous  eigenvalue  problem  (5.28),  which  have  been 
chosen  for  the  system  reduction.  The  calculation  of  these  de¬ 
rivatives  requires  a  high  numerical  effort  and  it  has  to  be  in¬ 
vestigated  if  their  influences  can  be  neglected.  It  can  be 
shown,  that  if  the  transformation  T  contains  all  possible  eigen¬ 
vectors,  that  means  no  system  reduction  is  achieved  (m  =  n), 
titese  terms  vanishes  exactly. 


IS  Noml 

nodes 

_ J 

20  Noml  nodes 

»g 

V,  -  426.5 

V,  -  419.6  n/s 

»*. 

analyt. 

nui. 

analyt. 

nua. 

1 

-6.084-3 

-7.350-3 

-8.117-3 

-8.668-3 

2 

7.971-4 

-6.534-4 

-5.751-4 

-2.127-3 

3 

1.593-4 

+2.742-4 

-2.117-4 

-4.145-4 

4 

1.457-1 

1.677-1 

1.602-1 

1.637-1 

Fig.  5.5:  Comparison  of  analytical/numeri^ 
flutter  gradients 
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In  Fig.  5.5,  the  results  for  the  derivation  of  the  flutter 
constraints  for  the  fin  -  example  described  in  [27]  are  de¬ 
picted.  The  comparison  between  the  numerical  and  the  ana¬ 
lytical  flutter  gradients  shows,  that  with  an  increasing  number 
of  modes,  the  quality  of  the  analytical  derivation  becomes 
much  better. 

6.  THE  SOFTWARE  SYSTEM  MBB- 
LAGRANGE 


Corresponding  to  the  Three-Columns-Concept  described 
in  section  2,  the  software  system  MBB-Lagrange  is  divided 
into  the  main  modules  structural  and  sensitivity  analysis,  op¬ 
timization  modeland  optimization  algorithms  (see  Fig.  6.1). 
The  proram  contains  design  models  for  cross-sections  or 
wall-thicknesses  of  isotropic  elements  for  layer  thicknesses 
and  layer  angles  of  fiber  composite  materials,  and  for  concen¬ 
trated  masses  as  well. 


Fig.  6.1:  General  program  architecture 


trol  pan,  a  pan  for  the  description  of  the  FE-model  and  one 
describing  the  optimization  model.  The  FE-description  is 
done  in  form  of  a  NASTRAN-Bulk-Data-Deck.  For  an  aero- 
elastic  analysis  the  aerodynamic  influence  matrices  must  be 
supplied  additionally. 

The  optimisation  results  are  documented  by  the  following 
data: 

•  a  file  for  describing  the  optimization  history, 

•  plotfiles  for  the  graphic  illustration  of  the  optimisation 
history, 

•  re-  and  warmstan  files  for  conanuing  an  optimisation, 

•  a  NASTRAN-Bulk-Data-Deck  of  the  optimised  struc¬ 
ture, 

•  an  IDEAS-Universal  file  of  initial  and  final  design  for  the 
graphic  illustration  of  the  strucmral  parameters,  displace¬ 
ments,  stresses,  strains,  values  of  the  safety  factor  etc. 

Besides  the  IDEAS-interface  there  is  also  an  interface  to 
the  pic-  and  post-processing-system  PATRAN. 

Since  many  different  optimisation  routines  are  available, 
a  user  must  either  define  them  "by  hand"  or  he  requires  a  se¬ 
lection  made  by  the  system.  In  this  case  a  rule-based  subproc¬ 
ess  [28]  will  send  some  questions  to  the  terminal  and  depend¬ 
ing  on  the  answers  and  the  information  on  the  design 
available  so  far.  a  heuristic  proposal  is  made.  A  user  may  ac¬ 
cept  the  proposed  method  and  parameters  or  he  may  choose 
another  code.  The  following  table  shows  an  example  for  the 
so-called  safety  factors,  which  indicates  if  a  strategy  will  be 
more  or  less  successful:  0  means  not  possible,  1(X)  means  it  is 
the  best. 

IBF  MOM  SLP  SRM  RQPl  RQP2  GRG  CONLIN 

0  63  64  0  70  70  14  75 


In  order  to  choose  the  most  suitable  optimization  algo¬ 
rithm  for  a  specific  problem,  the  following  algorithms  are 
supplied: 


1.  IBF 

2.  MOM 

3.  SLP 

4.  SRM 

5.  RQPI 

6.  RQP2 

7.  GRG 

8.  CONLIN 

9.  QPRLT 


Inverse  Barrier  Function, 

Method  of  Multipliers, 

Sequential  Linear  Programming, 
Stress  Ratio  Method, 

Recursive  Quadratic  Programming 
(Schittkowski), 

Recursive  Quadratic  Programming 
(Powell), 

Generalized  Reduced  Gradient, 
Convex  Linearization, 

Quadratic  Programming  with  Reduced 
Line  Search  Technique  (SQP-GRG- 
Combination) 


The  structural  and  sensitivity  analysis  consists  of  the  pro¬ 
cedures  for  determining  the  various  state  variables  and  their 
gradients,  to  characterize  the  static,  dynamic,  aeroelastic  and 
stability  behaviour. 

The  input  data  for  the  optimisation  are  divided  into  a  con- 


If  some  results  are  available  obtained  from  a  previous  run 
with  the  same  algorithm,  it  is  possible  to  perform  a  warm 
start,  i.e.,  continuation  of  the  iteration  which  was  interrupted 
before  by  exceeding  the  maximum  number  of  iterations.  Oth¬ 
erwise  a  cold  Stan  may  be  activated  starting  from  the  last 
computed  iterate  or  alternatively,  a  new  optimisation  cycle  is 
initiated  starting  from  the  originally  given  design  variables. 

MBB-Lagrange  possesses  a  very  flexible  failure  sytem 
and  it  is  out  of  the  scope  of  this  report  to  explain  all  of  its  fea¬ 
tures.  Severe  failures  interrupting  the  optimisation,  are  written 
to  a  output  file  and  are  sent  to  the  terminal.  By  activating  the 
failure  analysis,  a  user  will  see  the  same  failure  information 
again.  Subsequently  a  rule-based,  heuristic  proposal  of  a  suit¬ 
able  remedy  is  displayed  and  the  user  may  accept  the  pro¬ 
posed  action  or  not. 

7.  EXAMPLES 


The  following  sample  of  examples  gives  a  good  idea  of 
the  capabilities  of  MBB-Lagrange. 

Airbus  A300/600  Support  Beam 


T 


« 


Beside  the  stringers  and  flanges,  the  cabin  floor  influ¬ 
ences  the  mechanical  behaviour  of  the  fuselage  of  a  passenger 
or  cargo  aircraft  Important  parts  of  the  cabin  floor  are  a  large 
number  of  support  beams.  The  beams  are  connected  to  the  fu¬ 
selage  at  both  ends  and  are  supported  by  struts. 

The  suppon  beam  presented  here  has  a  shape  of  an  U- 
profile.  Because  of  the  symmetry  only  half  of  the  structure  is 
considered  for  the  FE-model,  which  gives  a  total  size  of  1068 
elements,  6409  degrees  of  freedom  and  one  load  case  (Fig. 
7.1). 


I 

I 


Fig.  7.1:  A300/600  Support  Beam 
(Finite  Element  Model) 

The  beam  is  manufactured  by  milling,  which  allows  a 
very  fine  discretisation  of  the  wall  thickness.  For  that  reason 
221  design  variables  could  be  defined.  The  stress  require¬ 
ments  are  assured  by  von  Mises  constraints  on  each  element. 

Additional  special  empirical  functions  for  compression 
(crippling  allowables),  according  to  the  german  aircraft  indus¬ 
try’s  design  manual,  are  defined  for  the  flanges  and  stiffening 
holes,  to  include  stability  constraints. 

The  optimisation  is  carried  out  by  the  SLP-algorithm  in  5 
to  13  iteration  steps  and  the  weight  is  reduced  by  about  30 
percent  depending  on  the  loading  condition. 

Fmine  of  a  Combat  Aircraft  Fuselage 

This  frame  is  located  in  the  inlet  for  the  engine  of  a  com¬ 
bat  aircraft.  It  is  a  typical  example  for  a  sizing  problem  of  a 
light  weight  structure  made  of  an  aluminium  alloy.  For  the 
formulation  of  the  optimisation  problem  it  is  important  to 
know  that  a  milling  machine  will  be  used  to  realize  a  variable 
thickness  distribution.  So  large  number  of  design  variables 
can  be  defined  in  order  to  calculate  the  optimal  thickness  of 
the  frame. 

The  finite  element  model  is  shown  in  Fig.  7.2.  It  involves 
975  degrees  of  freedom,  930  elements  and  97  load  cases. 


Fig.  7.2:  Finite  Element  Model  of  Frame 
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The  objective  function  is  the  weight  The  constraints  en¬ 
sure  static  requirements.  For  each  element  the  feasibility  of  the 
stresses  have  to  be  ensured.  Stability  constraints  are  taken  into 
account  to  prevent  buckling  and  lateral  instability  of  webs. 
Since  these  constraints  have  to  be  satisfied  for  each  single  load 
case  we  get  a  very  large  number  of  inequality  constraints  (m^ 

=  92829).  Together  with  the  thickness  design  variables  (n  = 
187)  it  is  a  really  large  scale  optimisation  problem. 

The  initial  design  has  infeasible  buckling  loads  and 
stresses.  The  c^timisation  algorithm  SLP  needs  20  iterations  to 
achieve  convergence.  Fig.  7.3  shows  the  c^timisatoin  history 
of  the  buckling  constraints,  where  respectively  the  most  criti¬ 
cal  constraint  of  the  corresponding  itoation  is  taken.  The  opti¬ 
mal  design  fulfills  all  static  requirements  of  all  97  load  cases 
and  achieves  a  weight  reduction  of  about  25  percent 
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Fig.  73:  History  of  most  critical  buckling 
constraint 

Horiawm  StaWlBcr  of  a  HeUcoBtcr  171 

The  structure  consists  of  an  airfoil  section  like  an  airplane 
wing  and  endplates  which  act  as  vertical  stabilizers.  The  upper 
and  lower  panels  of  the  airfoil  section  are  sandwich  plates  with 
a  honeycomb  core  and  aramid  fiber  face  sheets.  The  spar  is  an 
l-shaped  bar  with  sraps  made  from  unidirectional  carbon  fiber 
reinforced  material  and  a  shear  web  which  is  a  honeycomb 
sandwich  with  (TFRP  face  sheets. 


Fig.  7.4:  Design  variables  of  the  spar 
and  the  endplate 
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The  endpUtes  are  sandwich  plates  of  constant  thickness  con¬ 
sisting  of  aiamid  fiber  reinforced  face  sheets  and  a  honey¬ 
comb  core.  They  are  fixed  to  the  airfoil  section  by  screws. 

Fig.  7.4  shows  the  design  variables  of  the  span  and  the 
endplate. 

Three  load  cases  define  the  loading  of  the  stabilizer.  The 
structure  shall  withstand  these  loadings  with  a  factor  of  safety 
larger  than  l.S.  Sandwich  wrinkling  has  to  be  considered  as 
well  as  a  composite  failure  criterion  such  as  that  from  Tsai- 
Wu.  Additionally  to  these  constraints,  a  lower  bound  for  the 
first  eigenfrequency  is  given.  The  results  are  shown  in  Fig. 
7.5. 


•  10  12  14 
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Fig.  7.5:  Minimum  weight  of  the  horizontal 
stabilizer  versus  the  first 
eigenfrequency 

The  original  design  with  spar  has  a  structural  weight  of 
10.6  kg  (one  half  of  the  airfoil  section  plus  one  endplate)  and 
a  first  eigenfrequency  of  14  Hz.  Considering  only  strength  re¬ 
strictions.  the  weight  can  be  reduced  by  nearly  3  kg,  but  in 
this  case  the  first  eigenfrequency  drops  signifrcamly.  If  the 
first  eigenfrequency  is  held  constant,  then  a  weight  reduction 
of  about  2  kg  is  possible. 

From  Fig.  7.4  it  can  be  clearly  seen  that  the  design  with  a 
spar  is  far  better  than  a  design  without  a  special  spar.  Con¬ 
cerning  the  thickness  of  the  airfoil  section  it  can  be  stated  that 
for  low  stiffness  (e.g.  a  low  first  eigenfrequency)  there  is 
nearly  no  difference  in  the  weight  between  a  thin  (12%)  and  a 
thicker  (15%)  profile.  Only  if  a  high  stiffness  is  required  (a 
high  first  eigenfrequency)  a  thicker  aerodynamic  profile  is 
useful. 

Comnoaitc  Fin 

Fig.  7.6  shows  the  structural  model  of  the  well  known 
MBB-Fin.  The  cover  skins  of  the  fin  are  made  of  cartKin  fiber 
laminate  with  four  different  fiber  orientations  in  the  stabilizer 
and  three  in  the  rudder.  The  inner  supporting  structure  is  real¬ 
ized  by  an  aluminium  honeycomb  core.  The  fin  is  supported 
at  the  connection  points  to  the  fuselage  and  the  sdf^ss  of 
the  fuselage  is  modeled  with  a  general  stiffness  element  As 
stadc  load  cases  the  aerodynamic  forces  of  five  different 
flight  conditions  (different  sideslips  and  rudder  deflections; 
subsonic  and  supersonic)  are  chosen. 


lUbiliscr 


eoBftcctioB  poistt  to  tbe  fapeUc* 

Fig.  7.6:  Structural  model  of  the  fin 


For  this  concept  study  1862  constraints  were  defined: 

•  Stress  limitation  (isotT(q>ic  elements  1 19/load  case) 

•  Limitation  of  failure  safety  (FRC.  252/load  case) 

•  Aeroelastic  efficiencies  -  Fui  (0.8) 

(Ma  =  1 .8  (750  kts)  -  Rudder  (0.5) 

•  Rutter  speed  -  530  m/s 

The  problem  consists  of  102  sizing-design  variables  (one 
independent  design  variable  for  every  layer  in  every  element). 

The  sizing  optimum  results  in  a  weight  of  42.3  kg  (= 
100%)  for  the  variable  sldn  weight 

By  introducing  the  layer  angles  as  additional  design  vari¬ 
ables  it  is  possible  to  define  a  lot  of  other  optimisation  mod¬ 
els. 

A  model  with  sizing  plus  7  layer  angle  variables  (one  de¬ 
sign  variable  is  assigned  to  each  layer  of  the  stabilizer  and  of 
the  rudder)  leads  to  an  optimal  weight  of  34.6  kg.  An  optimal 
weiht  of  25.3  kg  is  achieved  for  an  opumisation  model  with 
additional  84  layer  angle  variables.  Fig.  7.7  (one  design  vari¬ 
able  fix’  every  angle  in  each  element  and  a  linking  of  die  first 
and  the  third  layer).  This  weight  is  the  theoretical  lower  limit 
and  it  will  be  not  manufacturable  but  it  shows  the  high  poten¬ 
tial  on  weight  saving  possibilities  including  fibre  orientation 
as  design  variables.  But  it  is  also  obvious,  that  manufacturing 
requirements  has  to  be  considered. 


Fi(.  7.7:  Optimal  thkkiiess  dbtributkm  for 
the  design  model  witii  84  layer  angle 
variables 
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For  this  reason  a  development  was  initiated  which  in¬ 
cludes  these  manufacturing  informations  as  additional  con¬ 
straints  into  the  optimisation  noodel.  By  this  it  will  be  possible 
in  the  near  future,  to  have  the  optimisation  as  the  "driving 
part"  in  the  cotnplete  composite  design  process. 

Integrated  Fin  Dwign 

This  example,  using  the  same  structural  and  aerodynamic 
noodel,  is  an  approach  to  an  integrated  design  analysis  with 
not  only  structural  sizing  variables  t  but  also  three  additional 
aerodynamic  design  variables: 

•  taper  ratio  X 

•  aspect  ratio  A 

•  surface  area  S 

The  interesting  response  quantity  of  this  study  is  the  unit 
side  load  p  as  basic  flight  mechanic  design  requirement  for  a 
vertical  fin.  It  depends  on  the  aerodynamic  derivative  on 
the  surface  area  S  and  the  aeroelastic  efficiency  fl.  The  state 
variable  equations  for  this  multidiscipline  problem  can  be  for¬ 
mulated  in  a  generalised  form  as: 

p  =  Cp  n  S  (Flight  mechanics) 

Cp  =  fy^  (X,  A)  (Aerodynamic) 

Ti  =  fj  (X,  A,  S,  t)  (Structure/Aeroelastics) 

The  internal  coupling  of  the  system  is  given  by  the  first 
equation.  The  system  sensitivity  equations  can  be  formulated, 
using  the  method  proposed  by  [13j.  The  partial  derivatives  of 
the  state  variables,  which  will  be  provicted  by  the  individual 
disciplines  are  on  the  right  hand  side  of  the  system  sensitivity 
equations. 
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Fig.  lA  System  sensitivity  equations 

The  derivatives  with  respect  to  the  aerodynamic  design 
variables  are  done  by  using  the  finite  difference  method  with 


The  results  of  this  design  study  are  shown  in  Fig.  7.9). 
The  state  variable  p  is  plotted  for  all  finite  difference  sensitivi¬ 
ties  of  the  deagn  variable  X,  A,  S  and  for  the  optimised  ele¬ 
ment  thicknesses,  with  an  aeroelastic  efficiency  fin  require¬ 
ment  of  80  percent.  (The  stress  and  strain  constraints  coming 
from  fiv.'  static  load  cases  are  in  the  optimal  design  also  ful- 
filkd.)  The  best  integrated  design  solution  is  got  with  a  10  per¬ 
cent  reduction  of  aspect  ratio.  In  this  case  the  lateral  unit  load 
will  be  slightly  increased  and  the  weight  is  reduced  by  7.S  per¬ 
cent 


8.  CONCLUSION 


This  paper  presents  a  way  of  solving  design  tasks  in  the 
aircraft  development  process  using  structural  optimisation 
methods.  As  design  criteria,  requirements  on  the  static,  dy¬ 
namic  and  aeroelastic  behaviour  of  aircrafts  are  considered. 
The  analysis  procedures  for  the  various  state  variables  describ¬ 
ing  this  behaviour  ate  based  on  the  Finite-Element-Method. 
The  appUcation  of  this  program  dettxmstrates  that  the  design 
process  can  be  supported  very  efficiently  by  the  structural 
timisation  method.  Another  important  advantage  is  the  fact 
that  the  structural  optimisation  enables  to  achieve  technically 
optimal  design. 

In  order  to  optimise  real-life  structures  many  inqtortant 
procedures  and  methods  ate  combined  in  the  optimisation  sys¬ 
tem  MBB-Lagrange.  Many  further  developments,  however, 
roust  follow.  Since  an  optimal  design  has  to  fulfill  all  demands 
on  the  structure  simultaneously,  suitable  conopletions  and  ex¬ 
tensions  of  the  structural  and  sensitivity  methods  as  well  as  the 
optimisation  models  (local  and  global  stability,  heat  transfer, 
acoustics,  thermal  stresses,  flight  mechanics,  and  control, 
manufacturing)  are  furthermore  required. 

For  fiber  con^iosiie  materials  in  particular  the  characteris¬ 
tic  possibilities  and  requirements  of  manuftctnring  must  be  in¬ 
cluded  in  the  optimisation  process  in  order  to  guarantee  that 
optimal  designs  can  be  produced  efficiendy  by  fully  utilising 
die  design  potential. 
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STRUCTURAL  OPTIMIZATION  OF  AIRCRAFT 
C.  CORNUAULT  -  C.  PETIAU 

DASSAULT-AVIATION,  78  Quai  Marcel  Dassault,  92214  SAINT-CLOUD 


ABSTRACT 

A  general  survey  of  Dassault  experience  and 
knowledge  on  Aircraft  Design  with  Optimization 
Methods  Is  depicted. 

This  survey  results  from  compiling  the 
developments  and  the  results  already  worked  out 
and  already  presented  In  several  papers  by 
C.  PETIAU  and  Al. 

Part  I  gives  a  detailed  description  of  the 
methodology.  The  special  features  of  optimization 
with  composite  materials  are  shown.  The 
organization  of  design  resulting  from  use  of 
optimization  techniques  Is  described  and 
techniques  neighbouring  optimization  as  model 
adjustment  are  reviewed,  as  well  as  further 
developments. 

Part  II  Illustrates  this  methodology  by  an  actual 
case  study  of  an  aircraft  design  by 
Oassault-Avlatlon  with  relevant  examples  of 
structural  and  aeroelastic  optimization  on  carbon 
structures  of  a  wing  and  a  fin. 

PART  I  -  METHODOLOGY  OF  STRUCTURAL  OPTIMIZATION 

1  -  INTRODUCTION 

The  structural  optimization  technique  has  been 
a  routine  process  at  Dassault  since  the  late 
1970s.  It  has  been  applied  for  all  projects  from 
the  Mirage  2000  to  the  Rafale. 

In  the  past,  the  design  of  a  structure  was 
achieved  by  the  'fully  stressed  design'  process 
(FSO),  which  consists  of  Iterations  of  drawing  and 
analyses,  with  reinforcement  where  the  structure 
Is  not  sufficiently  strong  and  lightening  where 
there  are  strength  margins.  However,  where  the 
only  constraints  on  a  metallic  structure  are  those 
relating  to  strength  of  material.  It  has  been 
demonstrated  (see  Ref.  1)  that  this  approach  Is 
neither  optimal  (maximization  of  stresses  Is  not 
equivalent  to  weight  minimization)  nor  efficient 
for  the  design  process.  In  practise,  the  designer 
Is  completely  unable  to  predict  Intuitively  any 
solution  when  constraints  relating  to  flexibility 
(such  as  elgenfrequencles,  aerodlstorslon  and 
flutter)  or  to  the  ply  arrangement  of  composite 
materials  are  Involved. 

Therefore  we  consider  today  that  the  use  of 
mathematical  optimization  Is  compulsory  for  the 
design  of  aircraft. 

We  have  built  the  structural  optimization  tool 
around  the  Dassault  softwares  CATIA  and  ELFINI 
which  Include. 


(1)  the  well-known  CAD  tool  CATIA,  which  gives  us 
geometry  and  mesh  generation, 

(2)  static  finite  element  analysis  for  linear  and 
nonlinear  problems, 

(3)  static  aeroelasticity,  calculation  and 
management  of  loads, 

(4)  linear  dynamics  :  calculation  of  elgenmodes, 
harmonic  and  transient  responses, 

(5)  nonlinear  dynamics  :  Impact  and  crash 
analysis,  landing  gear  and  aircraft 
Interaction, 

(6)  unsteady  aeroelasticity,  flutter,  coupling 
with  flight  control  system, 

(7)  fatigue  and  crack  propagation  analyses, 

(8)  heat  transfer  and  thermo-elastic  coupling, 

(9)  acoustic  and  elastoacoustic  coupling. 

The  optimization  monitor  covers  most  of  these 
branches. 

The  system  works  on  request,  in  either  an 
interactive  or  a  batch  mode,  and  uses  a  common 
data  base  managed  automatically.  Some  of  the  main 
common  characteristics  of  the  branches  are  : 

-  topological  dialogue  for  mesh  and  all  data 
generation.  All  properties  as  connectors  between 
nodes  and  elements  (geometry  connection  with 
CATIA  surface  element  characteristics,  etc.)  are 
described  by  blocks  of  constant  properties  in  a 
space  of  Indices  referring  to  node  and  element. 
The  process  leads  to  very  clean  meshes  for  all 
types  of  structure  from  the  whole  aircraft 
meshes  to  tridimensional  analyses  of  fitting 
details, 

-  a  wide  range  of  possibilities  for  visualization 
of  Inputs  and  outputs,  many  of  'wire  frame'  and 
'pixel'  types  of  pictures  for  displacement 
stresses,  failure  criteria  and  for  optimization 
design  variables,  active  constraints  and  safety 
margin  plots, 

-  advanced  mathematical  solution  :  the  solution  of 
linear  problems  Is  run  by  a  very  powerful 
variant  of  the  Frontal  Gauss  method,  which 
minimizes  the  computer  time  for  classical  linear 
problems. 

For  large  three-dimensional  problems 
the  use  of  the  conjugate  gradient  technique 
enables  the  same  level  of  performance  to  be 
maintained,  taking  Into  account  the  contact 
non linearities. 

For  geometric  nonlinear  problems  (membrane 
effects,  post-buckling,  snap  through,  etc.)  an 
original  algorithm  called  'preconditioned  BFGS 
with  exact  line  search'  has  been  developed.  This 
algorithm  benefits  directly  from  the  biquadratic 
character  of  total  potential.  It  can  handle  the 
most  severe  snap-through  conditions 


3-2 


which  shows  calculation  of  post-buckling  of  a 
curved  stiffened  panel  In  carbon  epoxy  material). 

We  must  underline  the  strong  practical 
Interest  of  the  post-buckling  analysis,  which 
enables  the  design  of  thin  composite  skin,  which 
buckles  before  ultimate  loading. 

We  are  going  to  present  a  more  detailed  view 

of  : 

-  the  optimization  technique  which  Is  mainly  used 
to  set  the  general  dimensions  of  the  structure. 
It  is  supported  by  FE  models  of  the  whole 
aircraft,  which  are  elaborated  only  from  the 
rough  definition  of  external  shape  and  Internal 
architecture,  the  result  of  this  optimization 
being  the  starting  point  of  detail  drawing, 

-  the  checking  analysis  which  comes  with  detail 
drawings, 

-  the  organization  of  drawing  and  analysis  which 
are  a  necessity  of  composite  design,  and  are 
present  possibilities  of  computer  tools. 

2  -  THE  OPTIMIZATION  METHOD 

We  present  the  operational  tool  as  it  was  used 
for  the  design  of  the  Rafale,  the  organization  Is 
iterative,  and  the  flow-chart  Is  shown  In  Fig.  1  : 

2.1  -  Cost  function 

The  current  goal  In  optimization  is  weight 
minimization.  Nevertheless,  1n  some  rd'.es,  weight 
can  be  taken  as  a  constraint,  t^'  n  joctive  being 
maximization  of  the  safety  margin  see  table  1). 

2.2  -  Design  variables 

The  characterization  of  the  optimisation 
design  variables  1'  made  on  groups  of  finite 
elements  (FE)  .  The  choice  of  these 

variables  partly  takes  into  account  manufacturing 
constraints  avd  tooling  rules  for  metallic 
material . 

For  a  composite  material,  the  design  variables 
are  the  number  of  plies  In  each  direction  for  each 
group. 

The  number  of  design  variables  often  reaches 
500,  which  can  act  simultaneously  over  several 
analysis  models. 

r .3  -  Constraints 

Constraints  Inequalities  come  from  the 
different  analysis  branches  of  ELFINI.  We  can 
consider  simultaneously  : 

(1)  various  failure  criteria  (Including  composite 
materials),  computed  from  static  stresses  for 
all  the  dimensioning  cases  of  loads, 

(2)  local  buckling  criteria, 

(3)  limited  displacements, 

(4)  aeroelastic  variation  of  aerodynamic 
derivatives, 

(5)  dynamic  natural  frequencies, 

(6)  flutter  speed  and  aeroelastic  dynamic  damping, 

(7)  various  technological  constraints  (such  as 
minimum  values  of  design  variables,  and 
limitations  of  the  thickness  variation  between 
adjacent  design  variables). 


The  constraints  considered  during  the  same 
optimization  can  come  from  several  analysis  models 
(e.g.  symmetric  and  antisymmetric  FE  aircraft 
model,  local  buckling  analysis  by  the 
Rayleigh-Ritz  method,  local  refined  FE  analysis, 
different  external  store  configurations  for 
dynamics  and  flutter,  variation  of  shape  because 
of  control  surface  deflections,  etc.). 

2.4  -  Sensitivities 

We  define  'sensitivities'  as  the  derivatives 
of  constraints  in  the  function  of  design 
variables.  The  principle  of  ELFINI  optimization  Is 
to  compute  these  derivatives  by  a  correct 
mathematical  process.  It  can  easily  be 
demonstrated  (see  Table  2  and  Refs  1  and  2)  that 
the  computation  of  derivatives  of  static  stresses, 
displacements,  and  aeroelastic  coefficients  is 
equivalent  to  solutions  with  a  'dummy'  case  of 
loads. 

The  number  of  loads  in  this  dummy  case  is  : 

(a)  number  of  loading  cases  x  number  of  design 

variables  if  formula  (1)  of  Table  1  Is  used  : 

(b)  number  of  constraints  if  formula  (2)  is  used. 

For  practical  problems  the  number  of  loads  in 
the  dummy  case  currently  reaches  several 
thousands,  and  their  solution  makes  up  the  main 
part  of  the  computer  cost  of  optimization. 

For  nonlinear  constraints  relative  to  the 
static  displacements  (equivalent  stresses,  failure 
criteria...).  The  operator  [J)Q'/,iX]ls  linearized 
near  X. 

When  constraints  are  eigenvalue  or  are  direc¬ 
tly  related  to  eigenvalue  (E.g.  eigenfrequency, 
linear  buckling  load,  divergence  or  flutter  speed, 
aeroelasticity  damping)  the  cost  of  their  deri¬ 
vation  is  negligible  (see  Tables  3, 4, 5, 6  and  7) 
However,  we  must  underline  that  these  derivations 
need  a  far  more  accurate  calculation  of  eigen¬ 
vector  than  those  needed  for  eigenvalue  analysis. 
Also,  we  have  found  that  It  is  very  difficult  to 
compute  with  reasonable  accuracy  derivations  of 
solution  of  problems  treated  with  the  classical 
modal  basis  reduction  (e.g.  dynamic  response, 
aeroelasticity),  in  practise  it  would  be  necessary 
to  compute  the  correct  mathematical  derivative  of 
the  basic  vectors.  This  is  mainly  why  we  have 
developed  a  static  aeroelasticity  approach  without 
the  basic  truncation  effect  (see  Ref.  1),  as  it 
leads  to  a  mathematically  exact  and  low-cost 
calculation  of  derivatives. 


2.5  -  Mathematical  optimization 

Starting  of  the  analysis  and  derivation  of 
constraints,  we  use  an  explicit  nonlinear  approxi¬ 
mation  of  the  constraints  In  terms  of  the  design 
variables,  mainly  the  formulation  In  Inverse 
variables.  Taking  as  new  variables  the  Inverses  of 
design  variables,  leads  to  minimization  of 
homographic  function  (weight)  subject  to  linear 
Inequalities.  This  problem  Is  easily  solved  by  a 
projected  conjugate  gradient  algorithm  (see  TABLE  8  ) 
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Sub-Ueratlon  process 

The  right  results  and  the  good  convergence  of 
our  algorithm  in  static  optimization  are  mainly 
due  to  the  explicit  approximation  of  constraints 
In  1/  . 

But  for  other  types  of  constraints  as  natural 
frequencies,  flutter,  speed,  dynamic  responses, 
this  explicit  form  in  l/;^i  has  no  theoretical 
basis  and  on  some  cases  we  could  have  a  bad 
convergence. 

Since  the  cost  of  the  calculation  of  these 
dynamic  constraints  Is  raised  relatively  low 
compared  to  the  analysis  cost,  it's  interesting  to 
carry  out  a  sub-iteration  process  In  order  to 
Improve  convergence. 

The  convergence  of  the  sub-iteration  process 
In  ensured  by  move  limits  and  a  relaxation  on 
admissible  value  of  the  constraints  making 
possible  the  detection  of  unfeasible  approximate 
problems. 

The  cost  of  the  mathematical  optimization  step 
is  low.  The  mathematical  optimization  step  gives  a 
prediction  of  the  optimum,  from  which  we  start  new 
Iterations.  The  number  of  iterations  needed  to 
obtain  global  convergence  ranges  from  three  to 
five  .  The  cost  of  all  of  the 
Iterations  of  optimization  ranges  from  about  eight 
to  15  times  the  cost  of  the  analysis. 

2.6  -  Final  touches 

Generally,  the  theoretical  optimum  obtained 
from  the  optimization  algorithm  needs  some 
modification,  as  it  often  does  not  represent  a 
realistic  design.  Starting  from  the  table  of 
constraint  derivatives,  the  final  touches  consist 
in  examining  Interactively  the  effect  of  small 
modifications,  made  directly  by  the  designer 
during  the  drawing.  The  program  instantaneously 
shows  the  new  safety  margin  and  any  violated 
constraints  . 

We  can  also  interactively  rerun  the 
mathematical  optimization  step  after  changing  the 
assigned  values  of  constraints. 

3  -  SPECIAL  FEATURES  OF  OPTIMIZATION  WITH 

COMPOSITE  MATERIAL 

The  organization  described  above  is  well 
suited  for  a  composite  material,  with  the  addition 
of  the  following  specifities. 

3.1  -  Failure  criteria  analysis  and  derivation 

Inside  the  optimization  loop  we  use  failure 
criteria  of  the  'TsaY-HIll'  family  : 


c  =  , 

/(# 

+  sr 

^jrad  / 

where 

Oj,ay 

and 

are  stress 

tensor 

components, and  and  SI  -  0  or  1  are 
criteria  parameters. 

The  arguments  of  the  criteria  are  adapted  to 
each  situation  (e.g.  tension,  compression, 
bending,  holed  panel,  etc.),  by  calibration  with 
more  sophisticated  criteria  and  test  results. 


Because,  at  a  given  point,  the  final  failure 
mode  is  not  known  beforehand,  it  Is  necessary  to 
handle  constraints  on  all  potential  failure  modes 
simultaneously.  This  Is  achieved  at  a  relatively 
low  cost  if  the  derivation  is  performed  in  two 
steps  : 

(1)  the  strain  tensor  and  its  derivative  are 
computed  by  formula  Cl)  of  Table  2  (three 
components  common  to  all  plies  with  membrane 
assumption), 

(2)  starting  from  the  strain  tensor  and  Hooke's 
law  for  the  material,  the  failure  criteria  and 
their  derivatives  are  calculated  ply  by  ply. 

3.2  -  Local  buckling  criteria 

Even  if  optimization  can  handle  global 
buckling  directly,  for  management  and  cost- 
effectiveness  it  is  generally  preferable  to 
calculate  and  derive  local  buckling  criteria  with 
the  following  post-processing  analysis  : 

(1)  using  the  general  FE  model,  stress  flows  of 
structural  meshes  are  calculated  and  derived, 

(2)  local  buckling  load  factors  and  their 

derivatives  are  calculated  by  a  Rayleigh-Rltz 
method  (see  Table  A). 

Sizes  of  meshes  for  local  buckling  analyses 
are  Independent  of  their  representation  in  the 
global  FE  model,  and  they  can  be  tuned  to  suit  the 
actual  stiffening. 

In  the  optimization  loop,  stacking  sequences 
are  not  taken  into  account  (It  is  assumed  that  the 
material  is  homogeneous  through  the  panel 
thickness),  for  the  sake  of  algorithm  simplicity, 
and  because  of  difficulties  in  expressing  the 
drawing  constraints  due  to  restrictions  of  cutting 
and  stacking  the  fiber  layers. 

The  order  of  buckling  modes  can  change  between 
iterations,  this  can  cause  a  non-convergence  of 
iterations  if  all  potential  buckling  modes  are  not 
controlled  simultaneously  (see  Ref.  2). 

3.3  -  Design  constraints 

These  constraints  express  the  fact  that  the 
results  of  optimization  must  correspond  to  a  real 
drawing  of  a  composite  panel,  which  must  be  made 
of  stacked  layers  with  special  rules  for  easy 
manufacture.  Design  constraints  are  handled  at  two 
levels  : 

(1)  inside  the  optimization  loop,  by  placing 
constraints  on  the  minimum  number  or  a  given 
minimum  proportion  of  plies  in  each  direction, 
or  on  a  maximum  slope  of  thickness  (these 
constraints  correspond  to  linear  inequalities 
In  design  variables), 

(2)  after  mathematical  convergence,  automatic 

rounding  of  thicknesses  is  used  to  obtain  a 
whole  number  of  plies,  and  a  special 

half-interactive  program  transforms  the 

stacking  of  plies  by  area,  which  are  the  rough 
output  of  optimization,  into  a  proper  lay-up 
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4  -  MULTIMODEL  OPTIMIZATION 

Optimization  has  to  provide  the  single 
physical  characteristics  of  a  structure  and  must 
take  all  the  sizing  considerations  Into  account. 
So  many  F.E.  models  for  other  type)  are  necessary, 
depending  on  the  studied  phenomena.  So 
optimization  must  Insure  : 

-  Identification  between  design  variables  defined 
on  several  models, 

-  data  transfers  between  models  (characteristics, 
boundary  conditions,  loads...), 

-  management  of  calculation  and  derivation  of 
constraints  defined  on  several  models, 

-  the  linking  of  all  the  design  variables  and 
constraints  (values  and  derivatives)  In  the 
single  explicit  optimization  step  giving  the 
optimum. 

The  organization  (see  TABLE  9  )  of  the  nwdels 
has  needed  some  software  Investment  and  Is  able  to 
manage  several  FE  meshes  with  several  boundary 
conditions,  mass  configurations  (modal  and  flutter 
analysis),  Mach  number  (aeroelasticity  and  flutter 
analysis).  Other  models  are  used  for  panel 
buckling  analysis. 

5  -  CHECKING  ANALYSIS 

It  must  be  understood  that.  If  an  optimization 
tool  Is  essential  to  achieve  a  good  general 
drawing  rationally,  the  result  must  be  justified 
In  detail,  using  more  complex  analyses  than  those 
which  can  be  handled  Inside  the  optimization  loop. 
The  most  typical  of  these  checking  analyses  are 
the  following  : 

(1)  effect  of  local  loads  (e.g.  fuel  tank 
pressures,  vibration,  thermal  load,  etc.), 

(2)  local  fatigue  analysis. 

(3)  damage  tolerance  analysis, 

(4)  detailed  local  analysis  of  holed  composite 
panel  (e.g.  point  stress  analysis), 

(5)  post-buckling  analysis  (see  Refs  3  and  4). 

Design  constraints  corresponding  to  these 
detail -checking  analyses  have  been  simplified  to 
be  handled  by  general  optimization.  These 
simplified  assumptions  must  be  validated  by  local 
checking  analysis. 

Effects  of  calibration  of  these  constraints 
can  be  examined  with  a  Lagrange  multiplier  of 
active  constraints  (handled  Interactively  by 
■final  touches'  modules)  or  by  re-running  the 
mathematical  optimization  step. 

6  -  ORGANIZATION  OF  DESIGN  PROCESS 

We  now  have  the  following  organization  for  the 
design  of  composite  structures,  from  the 
preliminary  project  to  the  delivery  of 
manufacturing  drawings  : 

(1)  start  from  a  CATIA  drawing  of  the  external 
shape  only  and  a  brief  definition  of  the 
Internal  architecture, 

(2)  elaboration,  by  CATIA-MESH,  of  a  first  simple 
general  FE  mesh  of  the  whole  aircraft 
(10-30000  dof)  with  approximate  cross-sections 
and  thicknesses  (see  Fig.  1).  The  model  Is 
adjusted  with  simple  cases  of  load, 

(3)  analysis  of  static  aeroelasticity  and  loads, 
which  give  the  envelope  cases  of  loads  and 
show  the  latent  problems  of  aeroelasticity. 


(4)  examination  of  Internal  load  fields  and 
stresses  for  selection  of  'strength  of 
material”  constraints  In  the  optimization. 

(5)  computation  of  dynamic  modes  with  the  various 

external  store  configurations  :  flutter 

problem  recognition, 

(6)  first  run  of  optimization, 

(7)  drawings  of  the  structure  supported  by  : 

(a)  an  Interactive  test  of  authoritative 
modifications  of  optimization  results  to 
make  drawing  easier,  together  with  use  of 
the  'final  touches'  module  ; 

(b)  changes  and  additions  of  constraints, 

(c)  critical  examination  of  'cost  of  require¬ 
ments',  directly  obtained  from  'Lagrange 
multipliers'  of  optimization.  This  allows 
appreciation  of  the  real  Influence  of  the 
safety  margin  of  uncertain  criteria 
(composite  materials), 

(d)  Detail -checking  analyses  supported  by 
methods  described  above  in  Section  5. 
These  are  performed  taking  proper  boundary 
conditions  In  the  FE  model  for  the  whole 
aircraft  via  a  super-element  technique. 
Detail-checking  analyses  must  validate  the 
simplified  criteria  used  for  mathematical 
optimization,  otherwise,  optimization  must 
be  re-run  with  calibrated  criteria. 

Although  a  single  optimization  run  lasts  only 
one  night,  the  optimization  job  can  remain  inside 
the  computer  for  more  than  6  months,  for  examina¬ 
tion  of  the  detail  analysis  effects,  the  influence 
of  the  choice  of  constraints  and  alternative 
designs . 

7  -  TECHNIQUES  NEIGHBOURING  OPTIMIZATION,  IDENTI¬ 
FICATION  AND  COMPUTATION  WITH  UNCERTAIN  DATA 

The  solution  of  these  problems  can  be 
considered  because  of  the  possibilities  of 
derivative  elaboration. 


7.1  Model  adjustment 

Generally,  this  involves  adjustment  of  the  FE 
dynamic  model  to  measured  natural  modes,  the 
unknowns  are  design  variables  of  local  thickness 
and  mass,  modal  deformation  and  frequencies.  The 
modal  equation  appears  as  an  equality  constraint, 
and  the  objective  is  to  minimize  the  'distance' 
between  the  measured  and  the  computed  modes.  The 
method  does  not  require  knowledge  of  the 
connection  between  computed  and  measured  modes, 
some  results  of  this  technique  applied  to  the 
Mirage  III  NG  are  shown  in  Fig.  2. 

For  general  cases  of  model  adjustment,  we  use 
a  simpler  technique.  The  objective  function  is  to 
minimize  the  'distance'  between  a  design  variable 
and  its  theoretical  value,  we  take  as  constraints 
the  fact  that  the  computations  must  give  measure¬ 
ments  with  a  given  approximation  (which  can  be 
objectively  estimated  from  the  accuracy  of  measu¬ 
rements). 

The  advantage  of  this  technique  over  the 
classical  mean-square  method  Is  that  under- 
determination  Is  not  possible,  if  a  design 
variable,  or  a  combination  of  design  variables,  is 
not  'observable'  by  measurement  the  process  gives 
the  theoretical  values  automatically.  In  Ref.  5  a 
good  example  of  this  process  for  flight  Identi¬ 
fication  of  aerodynamic  loads  was  given. 
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7.2  -  Computation  with  uncertain  data 

Sometimes,  at  the  start  of  a  problem,  the  data 
are  imprecisely  known,  the  idea  of  computation 
with  uncertain  data  is  to  find  the  "worst"  point 
in  the  uncertain  design  variable  space.  The 
problem  is  solved  by  two  approaches  ; 

(1)  find  the  "worst"  possible  point  by  minimiza¬ 
tion  of  a  safety  margin  function  inside  the 
authorized  space  of  design  variable 
variations,  and 

C2)  if  there  exists  a  possibility  of  failure, 
compute  the  probability  of  failure,  starting 
from  the  probability  density  of  the  design 
variables. 

We  have  now  started  to  apply  these  ideas  on 
flutter  and  vibro-ac  istic  analysis  of  preliminary 
projects. 

8  -  FURTHER  LEVELS  OF  OPTIMIZATION 

The  general  tendency  is  to  introduce  progres¬ 
sively  all  the  "arguments"  of  structural  design  in 
the  optimization  loop.  The  next  steps  of 
development  are  as  follows. 

8.1  -  Optimization  with  "bending"  design  variables 

This  does  not  give  rise  to  any  theoretical 
difficulties,  the  relative  complication  comes  from 
the  nonlinear  dependence  of  stiffness,  neutral 
surface  and  constraints  on  design  variables,  which 
complicates  program  writing. 

8.2  -  Optimization  with  post-buckling  analysis 

This  is  one  of  the  most  important  needs  of  the 
present  operational  optimization.  The  difficulty 
is  avoided,  generally  by  an  empirical  adjustment 
of  the  load  level  of  linear  buckling,  and  the 
results  of  optimization  are  checked  by 
post-buckling  analysis. 

The  correct  solution  is  not  much  more 
intricate  than  that  of  the  bending  case,  it  can 
easily  be  demonstrated  that  the  derivation  cost  is 
almost  that  of  the  linear  problem  ("dummy"  cases 
of  load  at  the  final  equilibrium  state)  (see 
table  10). 

8.3  -  Shape  optimization 

This  is  needed  in  many  pratical  problems  of 
varying  difficulty  (shape  of  stiffeners,  pressu¬ 
rised  vessels,  fitting,  etc.).  The  main  difficulty 
is  to  express  design  variables  and  "topological" 
constraints.  For  such  problems,  many  workers  and 
ourselves  have  elaborated  specimen  programs  for 
scholastic  cases,  but  for  a  really  operational 
tool,  it  is  necessary  to  introduce  geometrical 
design  variables  and  the  associated  "topological" 
constraints  at  the  level  of  a  CAD  system,  which 
requires  considerable  Investment. 

8.4  -  Optimization  in  heat  transfer  problems 

One  of  the  necessities  of  the  Hermes  project 
has  been  to  achieve  the  same  level  of  sophisti¬ 
cation  for  thermal  analysis  as  for  structural 
analysis,  we  have  met  the  need  for  a  thermal 
optimization  tool.  The  general  arrangement  of 
thermal  optimization  is  the  same  as  in  structural 
optimization.  The  complications  are  In  the 
transient  and  highly  non-11near  character  of 
thermal  problems. 


Fortunately,  it  can  be  demonstrated  that 
temperature  derivation  needs  the  solution  of  the 
same  differential  linear  equation  system  for  all 
design  variables,  and,  as  it  )s  integrated  at  the 
same  .  time  as  the  analysis,  it  does  not  need 
additional  factorization.  The  cost  of  derivatives 
is  therefore  relatively  lower  than  that  of 
solution  of  the  static  elasticity  (  "oblem. 

We  have  developed  a  joint  heat  transfer 
identification  process  with  computation  with 
uncertain  data,  which  is  needed  particularly 
because  of  the  random  or  badly  known  character  of 
many  data. 

8.5  -  Multidisciplinary  interactions 

For  a  combat  aircraft,  the  idea  should  be  to 
optimize  simultaneously  the  structure,  dimensions 
of  control  surfaces,  actuators  and  hydraulic 
power,  the  parameters  of  the  electrical  flight 
control  system,  and  the  aerodynamic  shape. 

This  state  of  grace  has  not  yet  been  reached, 
the  tendency  is  to  apply  optimization  to  each 
discipline  and  to  proceed  in  relation  to  the  other 
matters  by  a  "fixed  point  method"  or  by  simpli¬ 
fication  of  the  interactions.  Thus,  starting  from 
Lagrange  multipliers  obtained  from  the  optimi¬ 
zation  of  each  discipline,  it  is  possible  to 
"condense"  their  interactions,  for  instance,  as 
far  as  structure  is  concerned,  we  can  easily  give 
the  weight  cost  of  requirements  of  other  disci¬ 
plines  (exchange  rate  between  structure  weight  and 
roll  speed,  profile  relative  thickness,  etc.). 

9  -  CONCLUSION 

The  tendency  to  include  increasingly  detailed 
analyses  in  the  mathematical  optimization  loop  is 
hindered  by  the  difficulties  of  the  task.  The  tool 
described  above  represents  achievement  of  the 
first  level  of  structural  optimization,  where 
geometry  is  given  and  mass  and  sti^^fness  matrices 
are  linear  functions  of  design  variables,  signi¬ 
ficant  progress  is  not  easy.  It  corresponds  to 
including  in  the  optimization. 

(1)  "bending"  design  variables, 

(2)  nonlinear  and  post-buckling  analysis,  rules  of 
effective  width, 

(3)  stacking  order  of  plies  and  constraints  due  to 
cutting  of  layers  of  composite  material, 

(4)  shape  optimization,  which  is  also  implicitly 
necessary  in  the  above  functionalities. 

Apart  from  their  theoretical  difficulty,  these 
developments  need  a  higher  level  of  integration  of 
FE  optimization  with  CAD,  in  particular,  the 
architecture  of  the  CAD  system  must  support  the 
description  of  design  variables  and  of  drawing 
constraints. 

Another  promising  field  of  research  is  to  use 
expert  systems  to  pilot  the  design,  this  seems  to 
be  one  way  to  manage  optimization  with 
discontinuous  evolution  of  design  variables.  At 
present,  we  have  started  the  development  of  this 
technique  at  the  level  of  size  check  of  carbon 
fiber  panels.  It  rests  on  a  knowledge  base 
composed  of  rules,  relating  to  technological 
constraints  and  calculation  methods. 
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Fig.  1  Optimization  method  flow-chart 


TABLE  1 


TABLE  2 

Derivation  of  FE  Static  Solution 


MAXIMIZATION  OF  MARGINS 


Minimization  of  the  weight  is  not  the  only  aim 
of  the  optimization  of  structures.  Indeed  in 
projects,  by  example,  a  global  mass  is  allocated 
to  a  part  of  a  structure  (fuselage,  wing,  fin...) 
and  it's  interesting  to  design  the  safest 
structure  for  a  limited  mass. 

We  have  to  maximize  the  margin  for  a  limited 

mass 

Margin  :  jj_~  _  .51  -  1 

with  Cij  value  of  the  j  constraint 

admissible  value 

The  algorithm  has  to  minimize  the  yuj  = 
with  a  constraint  on  the  mass  : 

To  minimize  the  /ij  we  define  yu.  =  Max  ( /ij  ) 
and  the  problem  has  the  following  form  : j 
Min  yu 

for  constraints  concerned  by  the  margin 
Cj|^  for  other  constraints 

Mass  constraint 

In  reciprocal  variables  =  we  have  : 


(o<i  -  o<  io)  «  /I  <Ij 
(«i  -  ®<to)  <■  '^k 


Minimize  /i 


doj 

+  2 

1 

^o(  i 

dok 

+  1 

5dk 

i 

1 

oi 

4  ^0 

« i 


finite  element  analysis 

Displacement  computation: 

X  =  lATj-'F 

Strength,  stress  computation: 


Optimization  constraint  derivation 
Displacement  derivation: 

AJT  =  -iKl-'llA/flAT-  AFi 
Strength,  stress  derivation: 

Ao  =  -[^||j|X|-M|AAr|A'-AF|  (1) 

Ao  =  -[l/fl-JgJjllAA'lA'- AF|  (2) 

(1)  number  of  resolutions  =  number  of  load  cases 

(2)  number  of  resolutions  =  number  of  constraint  operators 


TABLE  3 

Derivation  of  Eigenvalues 


0  4  /t  «  1 


Analysis: 


_  eigenmodes:  K, 

This  optimization  problem,  with  a  linear 
objective  function  and  nonlinear  constraints  eigenvalues:  cu, 

(mainly  mass  constraint)  is  solved  by  a  sequence  ,  j 

of  linear  programming  using  moove  limits.  11^1  ~ 1^11  =  0 

Sensivity  analysis  of  eigenvalues: 

^\V'^\K\-o,]\M]\V,\  =  0 

2V',\\K\-a>^\M\\&.V,  +  FJIlA/f)  -tu?|AW||  +  Acu?  FJ | Af i I  K  =  0 
^  V]\\esK\  -  a>ntsM\\K 

2<o.V\\M]V, 
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TABLE  4 

Local  Buckling  Analysis  by  Rayleigh-Riiz  Method 


Rayleigh-Ritz  model 


External  load  fluxes:  0  = 


Normal  deflection: 
w  =  'La„„x‘  x'  •  Llx.y) 


0,0 


0yO 


=  P0O 


0,j.o| 


Buckling  factors 

Buckling  initiation: 


Ut,  =  W, 


M',  =  bending  elastic  energy 

fk'j  =  membrane  work  of  external  loads 


Pt/j) 


H',  =  pUi 

minp  =  WJUii=  =  )dW^ldV- pdUilBV  =  0 

|a:  -  pC|  =  0 


Derivation  of  buckling  factors 

^  y-KM-y 
^  y'-Gw-y 

dk  "  y^  cm-y  y'  G(0)  y 


TABLE  5 

-  Derivacives  of  extrema  in  transient  response 

Mechanical  equation  in  F.E.  basis  ; 

[M]  X  +  [k]  X  =  F  (t)  Q 

is  solved  in  the  reduced  basis  X  =  [v]  x 

by  integration  of  [m]  x  +  [k]  i  =  f  (l) 

with  [m]  =  [n]  V  [k]  =  [k]  V  F  (t)  =  F  (t) 

The  dynamic  stress  O'  can  be  written  ; 

[5Ci/:)x]x  =  [:)a/ix][v]  i=[i(i/5x]x 

Derivation  of  equation  6  gives 

[m]  ax  +[k]  ax  =-  [am]  X  -  [ak]  X 

which  is  similar  to  equation  (T)  except  the 
excitation.  If  this  dummy  excitation  can  be 
expressed  in  the  reduced  basis,  we  obtain  : 

[m]  Ax  +  [k]  Ax  =  -  [Am]  X  -[Ak]x 

[Am]=Vt[AM]v  [Ak]  »  Vt  [ak]  V 

On  the  extremum 

A  Q  (t)  «  [iO/Ox]  Ax^ 

d  <I(X,l)-[Aa/AX]d\  +  Xiig/5jrS'dt 


TABLE  6 


STATIC  AEROELASTIC  COEFFICIENTS 
ANALYSIS  AND  DERIVATIVES 


1  -  Basic  equations  of  static  aeroelasticity  with 
Finite  elements 


-  Discrete  pressure  field 

Kp  =  [  iKp/  Aqr]  qr  +  [  5Kp/  Aqs]  qs 

<1_  =  dgid  aerodynamic  effect  (incidence, 
control  surface  setting,  etc.) 
q  =  deformation  of  lifting  surface  expressed  in 
a  monomial  base. 

-  Flight  equation 

F  C.D.G.  =  i  pV^  ([Cr]  qr  +  [Cs]  qs) 

[Cr],  [Cs]  torseur  resulting  from  [AXp/Aqr], 
[  iKp/  iqs]  (aerodynamic  coefficients). 

-  Loads  on  F.E. 

[  :)F/  Aqr]  =  [R]  [  iXp/  iqr] 
F  =  I  P  V"  [R]  Kp  [  AF/  Aqs]  =  [R1  [  AKp/  Aqs] 

-F.E.  deformed  , 

,  [AX/  Aqr]  =  [K]-  [  AF/  Aqr] 

X  =  [K]-'  F  [  AX/  Aqs]  =  [K]-'  [AF/Aqs] 

-  Smoothed  transition  F.E.  qs  =  [L]X 

®  [A]l  =  [L]  [  AX/Aqr] 

qs  =  ipV'  ([Allqr  +  [A2]qs)  [A2]  =  [L]  [AX/Aqs] 

F.E.  =■  finite  elements 


Flexible  effect  elimination 


(D  qs  =  i  PV^Ju]  qr  ,  [pj  =  [D]'’  [Al]  , 

[D]  =  []  -  i  pr  [A2]] 


-  Flexible  aerodynamic  coefficients 
®  [c]  =  [Cr]  +  i  pv^  [Cs]  [jj] 

2  -  Derivative  relative  to  structural  parameters 

Differentiation  of  1  knowing  that 
A  ([K]'^  F)  =  -  [K]"''  [AK]  X 

Aqs  =  -  ^  PV"  [L]  [K]’’  lAK]  ([  A  X/  A  qr]  qr 
t  [  AX/  A  qs]  qs)  +  i  p  t  [A2]  Aqs 


By  eliminating  q^  from  2 

Aqs  =  -  J  p  V"  [D]’’  [L]  [K]’^  [AK]  C[  A  X/  A  qr] 

+  [  A  X/  Aqs]  [>i])  qr  =  [A>j]  qr 

Differentiation  of  3  [AC]  =  i  pv‘  [Cs]  [A)j] 

[AC]  =  i  PV^  [Cs]  [D]’’  [L]  [K]'’  [-AK] 

[[  AX/  Aqr]  +  [  AX/  Aqs]  [;j]] 

Preferable  in  the  form 

[AC]  =  i  PV"  CtK]"'  ([Cs]  [D]-'  [L])t)t[-A  K] 

[[  AX/  Aqr]  +  [  AX/  Aqs]  [>j]] 

which  means  solving  the  equilibrium  equation  for  a 
single  load  for  each  aerodynamic  coefficient  at  a 
given  Mach  and  dynamic  pressure 

[AC]  =  i  PV*  [Cs]  [D]’^  ([Xr^  [L]tV  [-AK] 

[[  AX/  Aqr]  [  AX/  Aqs]  [>i]] 

with  q  solutions  for  all  Mach  nuiftbers  and  dynamic 
pressures. 
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TABLE  7 

DERIVATION  OF  FLUTTER  SPEED  (J.P.  BREVAN'S  METHOD) 
1  -  Analysis 

-  We  atten^t  to  solve  the  flutter  equation  at  a 
given  Mach  number 

[K  (  \)  -  UJ^  (l+1g)*  M  (  X)  -  PV‘  A  q  =  0 


[K]  -  reduced  stiffness  matrix 
[M]  =  reduced  mass  matrix 
w  =  frequency  of  the  solution 
g  =  damping  fg  =  0  —e» flutter) 

V  =  velocity 

[A]  =  matrix  of  aerodynamic  forces 
P  =  left  solution 
q  =  right  solution 

-  In  the  simplified  form 

©  [D(  X  ,  w  ,  g,  V)]  q  =  0 

©  Pt  [DC  X  ,  g,  V)]  =  0 

with  g  =  0  Cor  given). 

2  -  Differentiation 

©  AC[D]q)  =  [AD]q  [D]  Aq  =  0 

by  multlplmylng  with  Pt 

_  Pt  [AD]q  PjtOAq  =  0 

we  obtain  ©  Pt  [AD]  q  =  0 

If  we  fix  the  damping  g 

0  Pt  [  0  D/  0X]q  ♦  Pt  [  :>  0/  Ou.  ]q 
+  Pt  [  0  D/  0  V]q  =  0 

This  complex  equation  gives  the  derivatives  of  the 
frequency  of  the  solution  and  the  flutter  speed. 


TABLE  8 

EXPLICIT  OPTIMIZATION 

The  main  Idea  Is  to  replace  the  exact 
formulation  fl{X)  of  the  constraints,  which  Is 
only  Implicitly  known  by  a  F.E.  analysis,  by  an 
explicit  approximation  C  * (X)  . 

<i**(X)is  selected  so  that  : 

-  at  X  =  X  analysis  Q  *C  X )  =  Q  E.F.  and 

id  "CX)/  iX  =  aa  E.F./6X 

-  Q  X )  must  be  exact  for  statically 

determinate  structure 

-  d  »C  X )  must  have  a  good  form  when  X  — »  0  and 
X  »  oo 

At  the  moment  we  think  that  the  most  efficient 
explicit  approximation  Is  In  reciprocal  variables. 

with 

•xl  =  1/Xi  alj  =  -l/<x1"  idj  E.F./:)X1 
d  i  "Cot)  =  aoj  +  1  alj  «xi 

aoj  =  d  E.F.  +  1  ai  j  C  «  1  -  cx  io) 
1 

The  optimization  problem  becomes  : 

1  ml /ot  1  minimal 

a^  j  +  I  ai  j  o<  1  0  j  ad  C «  ) 

«  4  1  /  Xi  mini 

With  constraints  on  local  buckling  criteria, 
the  admissible  values  are  not  constant  but 
function  of  parameters. 

The  explicit  optimization  problem  Is  solved 
with  a  conjuguate  projected  gradient  method 
Improved  by  an  efficient  normalisation  of  the 
tangent  Hesslen. 


Subject  to 
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TABLE  9 


TABLE  10 


MULTI -MODEL  OPTIMISATION 


OPTIMIZATION  OF  NONLINEAR  STRUCTURES 


Post-buckling  analysis  of  composite  structures 
is  one  of  the  most  significant  advance  of  the  last 
years  (Ref.  2).  The  next  challenge  is  the 
optimization  of  non-linear  structure  including 
post-buckling  behavior. 


Most  of  the  problems  can  be  solved  with  a 
similar  algorithm  to  those  of  linear  structures 
with  a  sequence  of  analysis  and  partial 
derivatives.  Derivatives  computation  are  relati¬ 
vely  less  expensive  than  in  linear  structures  : 


Pint  -  Fext  =0  Q  =  L.X 


dFint 

~rr 


Tx  =  - 


.liFint 
-1 


.  X  +  Ktg.  -^  =  0 

0  X 

^Fint  )(  5Q  _  . 

■  ’  5x  ■ 


( inm.  =  Ktg. 

'OX 
0  X 

ox 


But  this  type  of  algorithm  could  lead  to  bad 
convergence  on  post-buckled  structures  with 
snap-through  behavior. 


So  a  simultaneous  solution  of  the  analysis  and 
optimization  problems  can  be  considered. 


Min.  M  (  X ) 
a  (X,  X  )<  a  ad 

’^Rad  ’^tot  ° 


X  =  design  variables 

X  =  displacements 

^  =non-l inear  potential 


Recent  advances  in  minimization  method  based 
upon  preconditionned  matrices  and  explicit 
line-search  (Ref.  2)  will  be  intensively  used. 


INITIAL  MODEL 
F  •  8.52  Hz 


MESUREO  MODE 
WING  BENDING  2  NODES 
F  •  10,59  Hz 
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PART  II  -  ACTUAL  CASE  STUDY  OF  AN  AIRCRAFT 
DESIGNED  BY  DASSAULT  AVIATION 

EXAMPLE  OF  APPLICATION  OF  OPTIMIZATION  OF  CARBON 
EPOXY  STRUCTURE 

Ue  present  two  examples  of  optimization 
calculations  for  carbon  epoxy  parts  for  a  combat 
aircraft. 

1  -  Optimization  of  a  combat  aircraft  wing 

We  summarize  here  the  configuration  of  the 
optimization  of  a  carbon  epoxy  delta  wing  box, 
corresponding  to  the  mesh  presented  In  Fig.  1, 
with  the  design  variable  patch  of  Fig.  4  . 

We  have  used  two  analysis  models  for  static 
and  aeroelasticity  with  a  survey  of  flutter  on 
three  external  load  configurations  (see  Table  1). 

In  Fig.  5  we  present  the  history  of  conver¬ 
gence  In  weight.  Drawing  constraints  and  flutter 
constraints  have  been  successively  Introduced 
later,  to  study  their  Influence.  The  optimum 
values  of  design  variables  are  presented  In  Fig. 4 
and  the  corresponding  lay-up  of  plies  Is  shown  In 
Fig.  7  (obtained  automatically). 

In  Table  2  the  weight  sensitivities  of  wing 
panels  to  a  typical  project  hypothesis  obtained  by 
optimization  are  shown. 


TABLE  1 

Characterization  of  Optimization  of  a  Carbon  Epoxy  Wing 


Mo(M  / 

Model  2 

FE  models 

Wing  model  with  a 
representation  of  other  parts  of 
the  aircraft  by  super-element 
technique  (3544  doO:  symmetric 
and  antisymmetric  analysis 

Complete  plane  I3  003dor; 
symmetric  and 
antisymmetric  analysis 

Design  variable 

476  design  variables  (number  of  plies  in  four  directions) 

Static  cases  of 
loads 

24  cases  of  loads  combined 
from  symmetric  and 
antisymmetric  analysis 

0 

Failure  criteria 

476  failure  criteria  equivalent 
Tsaf'Hill  criteria 

Buckling  criteria 

144  critical  buckling  factors, 
obtained  from  77  local  buckling 
analyses  of  composite  plates  by 
Rayleigh- Ritz  method 

0 

Static  aeroelastic 
constraint 

0 

7  control  surface  efficiencies 
and  minimal  roll  speed 

Flutter 

5  flutter  speeds  and  60 
aeroelastic  dampings 
corresponding  to  three 
external  load  configuration 

Technological 

constraint 

374  constraints  on  composite  lay-up  (thickness  shape,  maximum 
and  minimum  ratio  between  each  ply  direction) 

TABLE  2 

influence  of  Design  Assumptions 

Design  hypothesis 

Weight 

(ratto) 

I  Composite  material  1 0 

Strength  of  material  constraints  only,  rough  from  computer 
optimization 


2  *  Aeroelasticity  constraint  I  19 

i  4  Aeroelasticity  4  technological  consuainis  I  23 

4  Weighifromnnaldetaileddraw)ng(reviewbycheckinganalyses)  136 

5  Aluminium  alloys  solution  2  10 

Strength  of  matenal  4  aeroelasticity  (comparable  with  line  3) 
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2  -  Optimum  design  of  a  vertical  fin 

The  layout  of  the  box  and  the  rudder  In  carbon 
epoxy  are  optimized  considering  the  static  load 
for  two  rudder  deflection  cases,  and  constraints 
on  rudder  aeroelastic  efficiencies  and  dynamic 
frequencies  (see  Fig.  6).  The  exact  configuration 
of  this  optimization  Is  shown  In  Table  3. 


TABLE  3 


Characteristics  of  Optimization  of  a  Vertical  Fin 


Model  1 

Model  2 

FE  models 

Rn  iitodd  (inOdoO  with  a  Fin  model  with  a  deflected 

super-element  of  the  whole  rudder 

aircrafi 

Design  variables 

237  design  variables  (number  of  fiber  layers  of  panels,  cross- 
sections  of  flanges  and  thicknesses  of  webs  for  spars  and  ribs) 

Stalk  load  cases 

3 

1 

Failure  criteria 

190  failure  criteria  on  composite 
maierUls  with  holes 

190 

Buckling  criteria 

91  buckling  criteria  computed 
from  Rayteifh-Riiz  models  for 
panel  buckling  analysis 

82 

Displacements 

1  displacement  on  the  step 
between  box  and  rudder 

0 

Aeroelasikity 

t  constraints  on  fin  and  rudder 
yaw  cfliciciicies  for  two  Mach 
numbers 

0 

Dynamk 

Frequencies  of  first  fiexion 
m<^  and  rudder  mode 

0 

Technology 

107  constraints  on  ply  distributions  and  on 

minimum  disance 

beffRKcn  lay-up  interruptions 
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Fig.  3  Post-buckling  analysis  of  a  curved  carbon  epoxy  panel 
(test  on  fuselage  panel  of  combat  aircraft) 


FI9.  4 


Optimization  of  a  carbon  epoxy  wing  :  (a)  upper  and  (b)  lower  panel  optimum  lay-up 


WEIGHT 

I 

«i 


2  3  4  5  6  7  8  9  10  11  12 


Iterations 


Iterations 


Iterations 

0  999  0  998  0  997  0998 

1024  0-998  0993  0999 


BUCKLING  CRITERIA 
(Raleigh  “Ritz  local  model) 


5  6  7  8 

iterations 


10  11  12 


Fig.  5  Optimization  of  carbon  epoxy  wing  :  history  of  convergence 
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Fig.  6  Optimization  of  a  carbon  epoxy  fin.  (a)  Model  1  :  no  rudder  deflection. 

(b)  Model  2  :  with  rudder  deflection,  (c)  Optimal  lay-up.  (d)  Active  constraints  at  optimum. 
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MULTIDISCIPLINARY  DESIGN  AND  OPTIMIZATION(^) 
by 

Jaroslaw  Sobieszczanski-Sobieski 
NASA  Langley  Research  Center,  MS  246 
Hampton  Virginia  2366S  U.S.A. 


SUMMARY 

Mutual  couplings  among  the  mathematical  models  of  physical 
phenomena  and  parts  of  a  system  such  as  an  aircraft  com¬ 
plicate  the  design  process  because  each  contemplated  design 
change  may  have  a  far  reaching  consequences  throughout  the 
system.  This  paper  outlines  techniques  for  computing  these 
influences  as  system  design  derivatives  useful  for  both  judg¬ 
mental  and  formal  optimization  purposes.  The  techniques  fa¬ 
cilitate  decomposition  of  the  design  process  into  smaller,  more 
manageable  tasks  and  they  form  a  methodology  that  can  easily 
fit  into  existing  engineering  organizations  and  incorporate  their 
design  tools. 

1.  INTRODUCTION 

The  engineering  design  process  is  a  two-sided  activity  as 
illustrated  in  Fig.  1.  It  has  a  qualitative  side  dominated  by 
human  inventiveness,  creativity,  and  intuition.  The  other  side 
is  quantitative,  concerned  with  generating  numerical  answers 
to  the  questions  that  arise  on  the  qualitative  side.  The  process 
goes  forward  by  a  continual  question-answer  iteration  between 
the  two  sides.  To  support  that  process  one  needs  a  compu¬ 
tational  infrastructure  capable  of  answering  the  above  ques¬ 
tions  expeditiously  and  accurately.  For  development  of  such 
an  infrastructure,  the  idea  of  “push  button  design”  ought  to 
be  discarded  in  favor  of  a  realistic  recognition  of  the  role  of 
human  mind  as  the  leading  force  in  the  design  process  and  of 
the  role  of  mathematics  and  computers  as  the  indispensable 
tools.  It  is  clear  that  while  conceiving  different  design  con¬ 
cepts  is  a  function  of  human  mind,  the  evaluation  and  choice 
among  competing,  discretely  different  concepts,  e.g.,  classical 
configuration  vs.  a  forward  swept  wing  and  a  canard  configu¬ 
ration,  requires  that  each  concept  be  optimized  to  reveal  its  full 
potential.  This  approach  is  consistent  with  the  creative  charac¬ 
teristics  of  the  human  mind  and  the  efficiency,  precision,  and 
infallible  memory  of  the  computer. 

The  computational  infrastructure  for  support  of  the  design 
process  entails  data  management,  graphics,  and  numerics.  The 
first  two  embodied  in  CAD/CAM  systems  ate  well-known  and 
are  taken  for  granted  as  a  framework  for  the  numerics.  The 
purpose  of  this  paper  is  to  introduce  some  imw  techniques 
which  may  be  regarded  as  a  subset  of  the  latter.  Included 
in  the  discussion  ate  the  system  behavior  derivatives  with 
respect  to  design  variables,  their  use  for  both  judgmental  and 
mathematical  optimization  purposes,  formal  decomposition 
of  a  system  imo  its  components,  and  ramifications  of  that 
decomposition  fw  system  sensitivity  analysis  and  optimization. 


('^Originally  pieaeiMed  under  the  title  “A  System  Approach 
to  Aircraft  Optimization"  as  Paper  No.  2  at  the  AGARD 
Workshop  on  Integrated  Design  Analysis  and  Optimization  of 
Aircraft  Structures,  1-2  May  1991,  Bath,  United  Kingdom. 


all  illustrated  by  aircraft  application  examples.  The  impact 
on  the  design  process  of  a  methodology  formed  by  these 
techniques  is  also  examined. 

2.  EFFECT  OF  DESIGN  VARUBLE  CHANGE  IN  A 
COMPLEX  SYSTEM 

An  aircraft  is  a  complex  system  of  interacting  parts  and  physi¬ 
cal  phenomena  whose  behavior  may  be  influenced  by  assigning 
values  to  the  design  variables.  Since  the  design  process  is,  gen¬ 
erally,  concerned  with  an  aircraft  that  does  not  yet  exist,  one 
works  with  its  surrogate — a  system  of  mathematical  models 
that  correspond,  roughly,  to  the  engineering  disciplines,  and  to 
physical  parts  of  the  vehicle.  These  mathematical  models  send 
data  to  each  other  as  depicted  in  the  center  of  Fig.  2,  and  they 
also  accept  design  variable  values  as  inputs  from  the  designers. 
To  know  how  to  change  these  design  variables,  designers  must 
know  the  answers  to  “what  if'  questions,  such  as  “what  will 
be  the  effect  on  the  system  behavior  if  the  design  variables 
X,  Y,  Z  will  be  changed  to  -t-  AX,  Y  -I-  AY,  X  +  AZ?”, 
implied  by  the  loop  in  Fig.  2. 

An  example  of  a  hypersonic  aircraft  in  Fig.  3  illustrates 
how  difficult  it  may  be  to  answer  an  “what  if’  question  for 
even  a  single  variable  change  in  a  complex  system  in  which 
everything  influences  everything  else.  Consider  a  structural 
cross-sectional  thickness  t  in  the  forebody  of  a  hypersonic 
aircraft  shown  in  the  upper  half  of  Fig.  3  as  a  design  variable 
that  is  to  be  changed.  The  lower  half  of  the  figure  depicts  a 
compfex  chain  of  influences  triggered  by  the  change  of  t  and, 
ultimately,  affecting  the  vehicle  performance.  The  change  of  t 
influences  the  position  of  the  bow  shock  wave  relative  to  the 
inlet  in  two  ways:  through  the  nose  deflection,  and  through  the 
weight  and  the  center  of  gravity  position  both  of  which  affect 
the  trimmed  angle  of  attack.  The  shock  wave  position  relative 
to  the  inlet  is  a  strong  factor  in  the  propulsive  efficiency  of 
the  engine  that,  in  turn,  combines  with  the  weight  to  influence 
the  aircraft  performance.  Additional  influence  on  performance 
is  through  the  angle  of  attack  whose  change  alters  the  vehicle 
aerodynamic  lift  and  drag.  The  resultant  modifications  of  the 
performance  may  require  resizing  of  the  vehicle  which,  of 
course,  may  be  a  sufficient  reason  to  change  t  again,  and  so  on, 
until  the  iteration  represented  by  the  feedback  loop  in  Fig.  3 
converges. 

The  above  iteration  engages  a  number  of  mathematical  mod¬ 
els  such  as  structures,  aerodynamics,  propulsion,  and  vehicle 
performance.  For  the  purposes  of  this  discussion,  each  such 
model  may  be  regarded  as  a  Mack  box  converting  input  to 
output  and,  consistent  with  the  black  box  concept,  the  inner 
workings  of  the  model  will  be  left  outside  of  the  scope  of  the 
discussion.  While  it  may  not  be  too  difficult  to  evaluate  the 
input-on-output  effect  for  each  single  Mack  box  taken  sepa¬ 
rately,  evaluation  of  the  resultant  change  for  the  entire  system 


4-2 


of  such  black  boxes  may  be  exceedingly  difficult,  especially 
when  iterations  are  involved.  In  general,  the  resultant  may  be 
a  small  difference  of  large  numbers,  so  even  its  sign  may  be 
impossible  to  predict  without  a  precise  reanalysis  of  the  entire 
system. 

To  generalize  from  the  above  example,  let  X  and  Y  denote 
the  system  input  and  output,  respectively,  e.g.,  the  structural 
cover  thickness  t  and  a  measure  of  performance  such  as  the 
aircraft  range.  Then,  the  derivative  dYfdX  is  a  measure  of  the 
influence  of  X  on  y  and  its  value  answers  quantitatively  the 
associated  “what  if’  question.  Mote  precisely,  the  derivative 
value  informs  only  about  the  rate  of  change  of  Y  at  the  value  of 
X  for  which  the  derivative  was  obtained.  Determination  of  the 
increment  of  Y  for  a  given  finite  increment  of  X,  if  Y{X)  is 
nonlinear,  can  be  done  approximately  by  a  linear  extrapolation 

(1)  y„ew  =  yoid  +  ^^ 

Capability  to  extrapolate  as  above  for  many  different  X  and 
Y  variables,  enables  one  to  decide,  either  judgmentally  or 
by  means  of  an  optimization  program,  which  variables  X  to 
change  and  by  how  much,  in  order  to  improve  the  design 
in  some  way.  However,  that  capability  is  predicated  on 
availability  of  the  derivatives  dY/dX  termed  the  system  design 
derivatives  (SOD).  For  large  system  analysis,  especially  if  the 
analysis  is  iterative,  its  is  advantageous  to  avoid  the  brute  force 
method  of  finite  differencing  on  the  entire  system  analysis  in 
computation  of  these  derivatives. 

2.1  System  Design  Derivatives 

Remembering  that  the  mathematical  model  of  an  engineering 
system  may  be  an  assemblage  of  a  large  number  of  mathe¬ 
matical  models  representing  its  components  and  the  governing 
physical  phenomena,  it  is  convenient  to  limit  the  discussion  to 
three  such  black  box  models  since  that  number  is  small  enough 
to  foster  comprehension  and,  yet,  large  enough  to  develop  a 
general  solution  pattern.  Ascribing  a  vector  function  repre¬ 
sentation  to  each  black  box,  the  set  of  equations  representing 
the  system  of  the  black  boxes  a,  0,  7  exchanging  data  as 
illustrated  in  Fig.  4  is 

Ya  =  Ya(X,Yff.Yy) 

(2)  Y0  =  Yg(X,Y„,Y^) 

Yj  =  Yy{X,Ya,Y0) 

The  y  and  X  variables  in  the  above  are  vectors  entered  in 
the  black  boxes  selectively,  e.g.,  some,  but  not  necessarily  all, 
elements  of  the  vectors  X  and  enter  the  black  box  0  as 
inputs.  Regarding  Y0{X,Ya,Yy)  as  an  example  of  a  black 
box,  the  arguments,  X,Ya,Yy,  ate  the  inputs  and  is  an 
output.  The  functions  in  eq.  2  are  coupled  by  their  outputs 
appearing  as  inputs,  hence  they  form  a  set  of  simultaneous 
equations  that  can  be  solved  for  Y  for  given  X.  The  act  of 
obtaining  such  a  solution  is  referred  to  as  the  system  analysis 
(SA).  In  the  presence  of  nonlinearities,  SA  is  usually  iterative. 

For  each  function  in  eq.  2,  one  can  calculate  derivatives  of 
output  with  respect  to  any  particular  input  variable,  assuming 


that  other  variables  are  fixed.  From  the  entire  system  perspec¬ 
tive,  these  derivatives  are  partial  derivatives  since  they  mea¬ 
sure  only  the  local  input-on-output  effect,  as  opposed  to  SDD 
which  are  total  derivatives  because  they  include  the  effect  of 
the  couplings.  To  prepare  for  furUier  discussion,  the  partial 
derivatives  corresponding  to  the  y-inputs  are  collected  in  the 
Jacobian  matrices  designated  by  a  pair  of  subscripts  identifying 
the  origins  of  the  output  and  input,  respectively.  For  example, 

(3)  m/dfYa] 

is  a  matrix  whose  y-th  column  is  made  of  the  partial  derivatives 
dYyifdYaj.  Assuming  the  length  otYfasN^  and  the  length 
of  yt  as  Na,  the  dimensions  of  matrix  J-ra  arc  Ny  X  Na¬ 
il  will  be  mnemonic  to  refer  to  the  partial  derivatives  in  the 
Jacobian  matrices  as  the  cross-derivatives. 

The  remaining  partial  derivatives  corresponding  to  the  X- 
inputs  are  collected  in  vectors,  one  vector  per  each  of  the  NX 
elements  of  the  vector  of  design  variables  X,  e.g, 

(4)  {dYa/dXt]'  =  [dYa/dXtl  k  =  l,...,NX; 

is  a  vector  of  the  length  Na  { '  denotes  transposition). 

Calculation  of  the  above  partial  derivatives  may  be  accom¬ 
plished  by  any  means  available  for  a  particular  black  box 
at  hand,  and  may  range  from  finite  differencing  to  quasi- 
analytical  methods  (ref.  1 ,  and  2). 

It  was  shown  in  ref.  3  that  differentiation  of  the  functions  in 
eq.  2  as  composite  functions  and  application  of  the  implicit 
function  theorem  leads  to  a  set  of  simultaneous,  linear,  alge¬ 
braic  equations,  referred  to  as  the  Global  Sensitivity  Equations 
(GSE),  in  which  the  above  partial  derivatives  appear  as  coeffi¬ 
cients  and  the  SDD  are  the  unknowns.  For  the  system  of  eq.  2, 
the  GSE  are 


•  I 

Jon 

(dYaldX^\ 

j'ay„/ax*'j 

I 

<  dYffjdX^ 

=  1  dYp/dXt 

I  . 

[  dYy/dXt 

[  dYy/dX„  J 

(5) 

These  equations  may  be  formed  only  after  the  SA  was  per¬ 
formed  for  a  particular  X,  a  particular  point  in  the  design 
space  because  the  computation  of  the  partial  derivatives  re¬ 
quires  that  all  the  X  and  Y  values  be  known.  For  a  given 
X,  the  matrix  of  coefficients  depends  only  on  the  system  cou¬ 
plings  and  is  not  affected  by  the  choice  of  X  for  the  right  hand 
side.  Hertce  that  matrix  may  be  factored  once  and  reused  in  a 
backsubstitution  operation  to  compute  as  many  sets  of  SDD’s 
as  many  different  Xi,  variables  are  represented  in  the  set  of 
multiple  right-hand-side  vectors. 

As  recommended  in  ref.  3,  numerical  solution  of  eq.  5  and 
interpretation  of  the  SDD  values  will  be  facilitated  by  normal¬ 
ization  of  the  coefficients  in  the  matrix  and  in  the  right  hand 
sides  by  the  values  of  and  Xo  of  the  y  and  JT  variables  for 
which  the  partial  derivatives  were  calculated.  The  normalized 
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coefficients  take  on  the  following  form,  illustrated  by  a  few 
examples  from  i-th  row  in  the  P  partition  in  eq.  5 

(6) 

where  the  normalization  coefficients  q  are 

^  -  ^ko 

Solution  of  the  normalized  eq.  S  yields  normalized  values  of 
the  SDD’s  from  which  the  unnormalized  values  may  always 
be  recovered  given  the  above  definitions. 

Formation  of  the  GSE  and  their  solution  for  a  set  of  SDD’s 
will  be  referred  to  as  the  System  Sensitivity  Analysis  (SSA). 

2.2  Utility  of  the  System  Design  Derivatives 
The  SDD  carry  the  trend  information  that  under  a  conventional 
approach  would  be  sought  by  resorting  to  statistical  data  or  to 
the  parametric  studies.  The  former  have  the  merit  of  capturing 
a  vast  precedent  knowledge  but  may  turn  out  to  be  ineffective 
if  the  vehicle  at  hand  is  advanced  far  beyond  the  existing 
experience.  The  latter  provide  an  insight  into  the  entire  interval 
of  interest  but  only  for  a  few  variables  at  a  time,  and  that  insight 
tends  to  be  quickly  lost  if  there  are  many  design  variables,  in 
which  case  the  computational  cost  of  the  parametric  studies 
also  may  become  an  impediment. 

In  contrast,  the  SDD  information  is  strictly  local  but  it  reflects 
the  influences  of  all  the  design  variables  on  all  aspects  of  the 
system  behavior.  Therefore,  the  SSA  should  not  be  regarded  as 
a  replacement  of  the  above  two  approaches  but  as  their  logical 
complement  whose  results  are  useful  in  at  least  two  ways. 

2.2.1  Ranking  design  variables  for  effectiveness 
A  full  set  of  SDD  for  a  system  with  NY  variables  in  Y  and 
NX  variables  in  A’  is  a  matrix  NY  x  NX.  The  j-th  column 
of  the  matrix  describes  the  degree  of  influence  of  variable  Xj 
on  the  behavior  variables  Y.  Conversely,  the  i-th  row  shows 
the  strength  of  influence  of  all  the  design  variables  X  on  the 
i-th  behavior  variable  Tj.  For  normalized  SDD’s.  comparison 
of  these  strengths  of  influence  becomes  meaningful  and  may 
be  used  to  rank  the  design  variables  by  the  ^gree  of  their 
influence  on  the  particular  behavior  variable.  This  ranking 
may  be  used  as  a  basis  for  judgmentally  changing  the  design 
variable  values  and  for  deciding  which  design  variables  to  use 
in  a  formal  optimization. 

An  example  of  such  ranking  is  illustrated  for  the  wing  of  a 
general  aviation  aircraft  shown  in  Fig.  5.  The  design  variables 
are  thicknesses  t  of  the  panels  in  the  upper  cover  of  the  wing 
box  and  the  behavior  variable  is  the  aircraft  range  R.  The 
chain  of  influences  leading  from  a  panel  thickness  to  the  range 
calculated  by  means  of  the  Breguet  formula  is  depicted  on 
the  left  side  in  Fig.  6.  In  the  Breguet  formula,  W,  denotes  the 
zero-fuel  weight  and  Wp  stands  for  the  fuel  weight.  Increasing 
t  in  one  of  the  panels  increases  the  weight  Wg  and.  in  general, 
reduces  the  drag  of  a  flexible  wing  by  stiffening  its  structure. 
Consequently,  the  range  is  influenced  in  conflicting  ways  that 
would  make  prediction  by  judgment  difficult.  However,  the 


corresponding  SSA  yields  the  SDD’s  for  the  upper  row  of  the 
wing  cover  panels  illustrated  by  the  heights  of  the  vertical  bars 
over  the  upper  wing  cover  panels  in  Fig.  6.  The  bars  show  that 
among  all  the  wing  cover  panels,  increasing  t  in  the  extreme 
outboard  panel  would  increase  range  the  most. 

2JL2  Gradient-guided  formal  optimization 
Most  of  the  formal  optimization  methods  applicable  in  large 
engineering  problems  use  the  first  derivative  information  to 
guide  the  search  for  a  better  design.  Since  the  SDD  values 
provide  such  information  for  all  the  Y  and  X  variables  of  in¬ 
terest,  the  SSA  may  be  incorporated,  together  with  SA,  in  a 
system  optimization  ptr>cedure  (SCH*)  based  on  the  well-known 
piecewise  approximate  analysis  approach  (e.g.,  ref.  4).  The 
SOP  flowchart  is  depicted  in  Rg.  7.  An  important  benefit  of 
the  SOP  organization  is  the  opportunity  for  parallel  processing 
seen  in  the  flowchart  operation  immediateiy  following  the  SA. 
In  that  operation,  one  computes  concurrently  the  partial  deriva¬ 
tives  of  input  with  respect  to  output  for  all  the  system  black 
boxes,  in  order  to  form  the  Jacobian  matrices  (eq.  3)  and  the 
right-hand-side  vectors  (eq.  4)  needed  to  form  the  GSE  (eq.  5) 
whose  solution  yields  the  SDD’s.  In  a  conventional  approach, 
these  SDD’s  would  be  computed  by  finite  differencing  on  SA. 
The  SDD  values  are  subsequently  used  in  Approximate  Anal¬ 
ysis  (extrapolation  formulas)  that  supplies  the  optimizer  (a  de¬ 
sign  space  search  algorithm)  with  information  on  the  system 
behavior  for  every  change  of  the  design  variables  generated  by 
that  optimizer,  and  does  it  at  a  cost  negligible  in  comparison 
with  the  cost  of  SA. 

A  generic  hypersonic  aircraft  similar  to  the  one  that  was  dis¬ 
cussed  in  Fig.  3  was  used  as  a  test  for  the  above  optimization. 
The  geometrical  design  variables  for  the  case  are  shown  in 
Fig.  8.  Additional  design  variables  were  the  deflections  of 
the  control  surfaces,  and  the  cross-sectional  structural  dimen¬ 
sions  of  the  forebody.  The  propulsive  efficiency  measured  by 
the  index,  defined  as  the  thrust  minus  drag  divided  by  the 
fuel  mass  flow  rate,  was  chosen  as  the  objective  function  to 
be  maximized.  The  aircraft  take-off  gross  weight  (TOGW) 
for  a  given  mission  is  very  sensitive  to  that  index,  thus  max¬ 
imization  of  the  iiKtex  effectively  minimizes  TOGW.  For  the 
reasons  discussed  in  conjunction  with  Fig.  3,  the  problem  re¬ 
quires  consideration  of  a  system  composed  of  aerodynamics, 
propulsion,  performance  analysis,  and  structures.  The  opti¬ 
mization  included  constraints  on  the  aircraft  as  a  whole  aixl  on 
behavior  in  the  above  disciplines.  Results  are  shown  in  Table  1 
in  terms  of  the  initial  and  final  values  of  the  design  variables 
(cross-sectional  dimensions  omitted)  and  of  the  objective  func¬ 
tion,  all  normalized  by  the  initial  values.  Considering  that  the 
initial  values  resulted  from  an  extensive  design  effort  using 
a  conventional  ^tproach,  the  nearly  13%  improvement  in  the 
propulsive  efficiency  was  regarded  as  very  significant  indeed. 

Another  example  of  the  SOP  application  is  the  case  of  a  hyper¬ 
sonic  interceptor  (Fig.  9a)  reported  in  ref.  5.  The  optimization 
objective  was  the  minimum  of  TOGW  for  the  mission  pro¬ 
file  illustrated  in  Fig.  9b.  The  system  comprised  the  modules 
of  the  configuration  geometry,  configuration  mass  properties, 
mission  performance  analysis,  aerodynamics,  and  propulsion 
as  depicted  in  Fig.  10,  and  the  design  variables  were  the  wing 
area,  scale  factor  for  the  turbojet  engine,  scale  factor  for  the 
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ramjet  engine,  and  the  fuselage  length.  The  constraint  list 
included  a  limit  on  the  time  needed  to  reach  the  combat  zone, 
the  take-off  velocity,  and  the  fuel  availaUe  mass  being  at 
least  equal  to  the  one  required  (the  fuel  balance  constraint). 
It  should  be  noted  that  in  a  conventional  approach  to  aircraft 
design,  satisfaction  of  the  latter  constraint  is  one  of  the  prin¬ 
cipal  goals  in  development  of  a  baseline  configuration  whose 
iminovement  is  subsequently  sought  by  parametric  studies  in 
which  the  design  variables  are  varied  while  always  striving  to 
hold  the  fuel  balance  constraint  satisfied.  In  contrast  to  that 
practice,  the  optimization  reported  in  ref.  S  allowed  the  fuel 
balance  constraint  to  be  violated  in  the  baseline  configuration 
and  achieved  satisfacticni  of  that  constraint  in  the  course  of  the 
optimization  process.  This  demonstrated  that  an  optimization 
procedure  may  do  more  than  Just  improve  on  an  initial,  feasible 
configuration;  it  can  actually  synthesize  an  optimal  configura¬ 
tion  starting  with  one  that  is  not  even  capable  of  performing  a 
required  mission. 

The  optimization  results  are  illustrated  by  a  vertical  bar  chart  in 
Fig.  I  i  that  shows  the  changes  of  the  design  variables  and  of  a 
significant  (13%)  improvement  of  the  objective  function.  The 
figure  shows  also  that  the  initially  violated  constraints  of  time 
to  intercept  and  take-off  velocity  were  brought  to  satisfaction 
in  the  optimal  configuration.  The  SOP  converged  in  only  4  to 
5  repetitions  of  SA  and  SSA. 

3.  MERITS  AND  DEMERITS 

Before  discussion  of  the  ramifications  of  the  above  sensitivity- 
based  optimization  in  a  system  design  process,  it  may  be  useful 
to  examine  briefly  the  merits  and  demerits  of  the  proposed 
approach  relative  to  the  conventional  technique  of  generating 
SSD  by  finite  differencing  on  the  entire  SA. 

3.1  Accuracy  and  Concurrent  Computing 
The  SSA  based  on  eq.  5  has  two  unique  advantages.  First, 
the  accuracy  of  SDD  is  intrinsically  superior  to  that  obtainable 
from  finite  differencing  whose  precision  depends  on  the  step 
length  in  a  manner  that  is  difficult  to  predict.  As  pointed  out 
in  ref.  6  it  is  particularly  true  in  the  case  of  an  iterative  SA 
whose  result  often  depends  on  an  arbitrary,  “practical”  con¬ 
vergence  criterion.  Second,  there  is  an  opportunity  for  con¬ 
current  computing  in  the  generation  of  the  partial  derivatives 
which  exploits  the  technology  of  parallel  processing  offered 
by  multiprocessor  computers  and  computer  networks.  Con¬ 
current  computing  also  enables  the  engineering  workload  to 
be  distrih-ited  among  the  specialty  groups  in  an  engineering 
organization  to  compress  the  project  execution  time. 

3J  Computatkmal  Cost 

Experience  indicates  that  in  large  engineering  applications, 
most  of  the  optimization  computational  cost  is  generated  by  the 
finite  difference  operations.  Therefore,  relative  reduction  of  the 
cost  of  these  operations  translates  into  nearly  the  same  relative 
reduction  of  the  cost  of  the  entire  optimization  procedure. 

The  computational  cost  of  the  SSA  based  on  eq.  3,  designated 
Cl,  may  be  reduced,  in  most  cases  very  decisively,  below 
that  of  finite  differencing  on  the  entire  SA,  denoted  by  Cj, 
but  to  achieve  that  reduction  the  analyst  should  be  aware  of 
the  principal  factors  involved.  To  define  these  factors,  let  the 


computational  cost  of  one  SA  be  denoted  by  CSA  while  CBAi 
will  stand  for  the  computational  cost  of  one  analysis  of  the 
t-th  black  box  in  the  system  composed  of  NB  black  boxes. 
The  t-th  black  box  receives  an  input  of  NXi  design  variables 
X,  and  NYi  variables  Y  from  the  other  black  boxes  in  the 
system.  Assuming  for  both  alternatives  the  simplest  one  step 
finite  difference  algorithm  that  requires  one  reference  analysis 
and  one  perturbed  analysis  for  each  input  variable,  the  costs 
Cl  aitd  Ci  may  be  estimated  as 

Cl  =  ^(1  +  NXi  +  NYi)CBAi-, 

(7)  < 

Cj  =  (1  -I-  NX)CSA 

Even  though  one  may  expect  CBAi  <  CSA,  a  sufficiently 
large  NYi  truty  generate  Ci  >  Cj  and  render  SSA  based  on 
eq.  5  unattractive  compared  to  finite  differencing  on  the  entire 
SA.  This  points  to  NYi,  termed  the  interaction  bandwidth, 
as  the  critical  factor  whose  magnitude  should  be  reduced  as 
much  as  possible.  Reducing  the  interaction  bandwidth  requires 
judgment  as  illustrated  by  an  example  of  an  elastic,  high  aspect 
ratio  wing  treated  as  a  system  whose  aeroelastic  behavior 
is  modeled  by  interaction  of  aerodynamics  and  structures, 
represented  by  an  CFD  analysis  and  Finite  Element  analysis 
codes,  respectively.  If  one  let  the  full  output  from  each  of 
these  black  boxes  be  transmitted  to  the  other,  there  might 
be  hundreds  of  pressure  coefficients  entering  the  structural 
analysis  and  thousands  of  deformations  sent  to  the  aerodynamic 
analysis.  With  the  NYi  values  in  the  hundreds  and  thousands, 
respectively,  it  would  be  quite  likely  that  Ci  >  C%.  However, 
one  may  condense  the  information  flowing  between  the  two 
black  boxes  by  taking  advantage  of  the  high  aspect  ratio  wing 
slenderness.  For  a  slender  wing  it  is  reasonable  to  represent 
the  entire  aerodynamic  load  by,  say.  a  set  of  S  concentrated 
forces  at  each  of  10  separate  chords,  and  to  reduce  the  elastic 
deformation  data  to,  say,  elastic  twist  angles  at  7  separate 
chords.  This  condensation  reduces  the  NYi  values  to  SO 
for  structures  and  7  for  aerodynamics.  In  the  finite  element 
code,  that  implies  SO  additional  loading  cases  all  of  which 
can  be  computed  very  efficiently  by  the  multiple  loading  case 
optiott — a  standard  feature  in  finite  element  codes.  The  CFD 
code  would  have  to  be  executed  only  7  additional  times.  Thus, 
the  advantage  of  the  interaction  bandwidth  condensation  is 
evident.  In  general,  a  condensation  such  as  the  one  described 
above  for  a  particular  example  may  be  accomplished  by 
the  reduced  basis  methods,  arttong  which  the  Ritz  functions 
approach  is,  perhaps,  the  best  known  one. 

3J  Potential  Singularity 

One  should  be  aware  when  using  SSA  based  on  eq.  S  that, 
in  some  cases,  the  matrix  of  coefficients  in  these  equations 
nuy  be  singular.  In  geometrical  terms,  a  solution  in  SA 
may  be  interpreted  geometrically  as  a  vertex  of  hyperplanes 
on  which  the  residuals  of  the  governing  equations  for  the 
black  boxes  involved  are  zero.  As  pointed  in  ref.  3,  eq.  S  are 
well-conditioned  if  these  hyperfflanes  inwrsect  at  large  angles, 
ideally  when  they  are  munially  orthogonal.  For  two  functions 
of  two  variables  the  zero-residual  hyperplanes  reduce  to  the 
zero-residual  contours,  and  an  example  of  a  nearly-orthogonal 
solution  intersection  is  shown  in  Fig.  12a.  In  some  cases. 
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the  intersection  angles  may  tend  to  be  very  acute,  in  the  limit 
they  may  be  zero  in  which  case  a  solution  exist  by  virtue  of 
tangency  of  two  curved  contours  as  illustrated  in  Fig.  12b.  It 
is  shown  in  ref.  3  that  eq.  S  imply  local  linearization  of  these 
comours  in  the  vidnity  of  the  intersection  point  so  that  the 
solution  point  is  interpreted  as  an  imersectim  of  the  tangents. 
Consequently,  in  the  situation  defncted  in  Fig.  12b  the  tangents 
coindde  and  the  matrix  of  eq.  S  becomes  singular.  In  such  a 
case,  eq.  S  should  be  replaced  by  an  alternative  formulation  of 
the  system  sensitivity  equations  in  ref.  3  based  on  residuals. 

There  were  no  cases  of  singularity  reported  so  far  in  any 
applications  probably  because  the  system  solutions  of  the  type 
illustrated  in  Fig.  12b  characterize  an  ill-posed  system  analysis 
usually  avoided  in  practice. 

3.4  Discrete  Variables 

Ndther  the  reference  technique  nor  the  SSA  based  on  eq.  5  can 
accommodate  truly  discrete  design  variables.  Truly  discrete 
design  variables  are  defined  for  the  purposes  of  this  discussion 
as  those  with  respect  to  which  SA  is  not  differentiable.  These 
are  distinct  from  quasi-discrete  variables  with  respect  to  which 
SA  is  differentiable  but  which  may  only  be  physically  realiz¬ 
able  in  a  set  of  discrete  values.  An  example  of  the  former  is 
an  engine  location  on  the  aircraft:  either  under  the  wing  or  at 
the  aft  end  of  the  fuselage.  An  example  of  the  latter  is  sheet 
metal  thickness  available  in  a  set  of  commercial  gages. 

In  the  case  of  truly  discrete  design  variables,  different  combina¬ 
tions  of  such  variables  define  different  design  concepts  (alter¬ 
natives)  and  each  concept  may  be  optimized  in  its  own  design 
space  of  the  remaining  continuous  variables,  in  order  to  bring 
it  up  to  its  true  potential.  Then,  one  may  choose  from  among 
the  optimal  alternatives.  Occasionally,  a  continuous  transfor¬ 
mation  might  be  possible  between  two  concepts  that  seem  to 
be  discretely  different.  For  example,  a  baseline  aircraft  with  a 
canard,  a  wing,  and  a  conventional  tail  may  be  reshaped  into 
any  configuration  feamring  all,  or  only  some  of  these  three 
lifting  surfaces.  This  is  so  because  a  sensitivity-guided  SOP 
may  eliminate  a  particular  feature,  if  a  design  variable  is  re¬ 
served  for  that  feature  and  if  the  feature  is  present  in  the  initial 
design  (however,  a  feature  initially  absent  cannot,  in  general, 
be  created). 

3.5  Noa-utUization  of  Dudpliiuury  Optimization 
Organization  of  the  SOP  discussed  above  may  be  described  as 
“decomposition  for  sensitivity  analysis  followed  by  optimiza- 
titm  of  the  emire,  undecomposed  system”.  It  may  be  regarded 
as  a  shortcoming  that  the  procedure  leaves  no  clear  place  for 
the  use  of  the  vast  expertize  of  optimization  available  in  the  in¬ 
dividual  Mack  boxes  representing  engineering  disciplines.  Ex¬ 
amples  of  such  local,  disciplinary  optimization  techniques  are 
the  t^mality  criteria  for  minimum  weight  in  structures,  and 
shaping  for  minimum  drag  for  a  constant  lift  in  aerodynamics. 
It  appears  that  combining  these  local,  disciplinary  optimiza¬ 
tion  techniques  with  the  overall  system  optimization  should 
benefit  the  latter.  Indeed,  one  way  in  which  theae  techniques 
may  be  used  without  changing  anything  in  the  SOP  oiganiza- 
tion  described  above  is  in  the  SOP  initialization.  Obviously, 
starting  SOP  ftum  a  baseline  system  composed  of  the  blacfc 
boxes  already  preoptimized  for  minimum  weight,  minimum 


drag,  maximum  propulsive  efficiency,  etc.  should  accelerate 
the  SOP  convergence  and  improve  the  end  result.  Such  local 
t^mizations  could  be  accomplished  separately  for  each  black 
box,  assuming  X  and  guessing  at  the  F  inputs. 

Beyond  that,  the  issue  of  incorporating  the  local,  disciplinary 
optimization  in  SOP  remains  to  be  a  challenge  for  further 
development.  Some  solutions  were  proposed  in  ref.  7  and  8  but 
their  effectiveness  is  yet  to  be  proven  in  practice.  This  issue 
will  be  taken  up  again  in  the  later  discussion  in  conjunction 
with  the  special  case  of  a  hierarcluc  system  decomposition 
which  does  accommodate  the  local  optimizations. 

4.  FORMAL  DECOMPOSITION 
When  the  system  at  band  contains  a  large  number  of  black 
boxes  and,  especially,  if  there  is  little  or  no  experience  with 
its  solution,  it  is  useftil  to  apply  a  formal  technique  to  deter¬ 
mine  the  data  flow  among  the  black  boxes.  The  data  flow 
information  is  useful  because  it  characterizes  the  system  as 
non-hierarchic,  hierarchic,  or  hybrid,  and  this,  in  turn,  helps  to 
choose  an  optimization  approach  and  to  establish  an  efficient 
organization  of  computing.  Such  formal  techniques  are  avail¬ 
able  in  Operations  Research  and  some  of  them  were  adapted 
for  the  system  analysis  and  optimization  purposes,  e.g.,  ref.  9. 

4.1  Af-square  Matrix 

A  brief  introduction  to  one  such  technique  begins  with  a 
formalization  of  a  black  box  (a  module)  in  the  system  as  one 
that  receives  inputs  through  the  top  and  bottom  horizontal  sides 
and  sends  the  output  through  the  left  and  right  vertical  sides  as 
as  shown  in  Fig.  13.  Using  that  formalism,  one  can  represent  a 
four-module  system  example  depicted  by  the  diagram  (known 
as  the  graph-theoretic  format)  in  Fig.  14a  in  a  different  format 
shown  in  Fig.  14b.  That  format  is  known  as  the  Af-$quare 
Matrix  format  because  N  modules  placed  along  the  di^onal 
fcrnn  an  table.  The  N-sqiiaie  Matrix  format  assumes  that 
the  modules  are  executed  in  order  from  upper  left  to  lower 
right  (although,  if  possible,  concurrent  executions  are  allowed). 
If  the  execution  order  is  not  yet  known,  the  order  along  the 
diagonal  may  be  arbitrary.  Referring  to  Fig.  13,  each  module 
may,  potentially,  send  data  horizontally,  left  and  right,  and 
receive  vertically  from  above  and  from  below.  The  actual 
data  transmissions  from  and  to  i-th  module  are  determined  by 
comparing  the  module  input  list  to  the  predecessor  module 
output  lists  while  moving  upward  in  column  t.  Wherever  a 
needed  input  item  is  found  on  the  output  list  from  module  j, 
a  dot  is  placed  at  the  imersection  of  the  i-th  column  and  y-th 
row  as  a  data  junction  indicating  transmission  of  output  from 
module  j  to  input  of  module  i.  After  the  predecessor  module 
search  gets  to  the  first  module,  it  switches  to  module  i  -H  1 
and  continues  downward  through  all  the  successor  modules  to 
module  N.  If  more  than  one  source  is  found  for  a  particular 
input  item,  a  unique,  single  source  must  be  judgmentally 
selected.  However,  an  output  item  may  be  used  by  several 
receiver  modules  and  nuiy  also  be  sent  to  the  outside.  The 
input  items  that  could  not  be  found  in  the  vertical  search  are 
designated  primaiy  inpuu  to  be  obtained  ftwn  the  outside  of 
the  system.  The  above  search  is  readily  implememable  on  a 
computer. 
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When  the  above  search  procedure  is  completed  for  all  the 
modules,  the  result  is  an  Af-square  Matrix  as  in  Fig.  14b  that 
conveys  the  same  information  as  the  diagram  in  Fig.  I4a  but 
is  amenable  to  computerized  manipulation.  To  see  what  such 
manipulation  may  achieve,  observe  that  each  dot  in  the  upper 
triangle  of  the  Af-square  Matrix  denotes  an  instance  of  the 
data  feedforward,  and  each  dot  in  the  lower  triangle  notes  an 
instance  of  the  data  feedback.  Of  course,  every  instance  of 
a  feedback  implies  an  iteration  loop  required  by  the  assumed 
diagonal  order  of  the  modules.  However,  that  order  may  be 
changed  at  will  by  a  code  that  may  be  instructed  to  switch 
the  modules  around,  with  the  associated  permutations  of  the 
rows  and  columns  to  preserve  the  data  junction  information,  in 
order  to  eliminate  as  many  instances  of  feedback  as  possible. 
If  all  of  them  are  eliminated  the  system  admits  a  sequential 
module  execution,  and  may  offer  opportunities  for  concurrent 
executions  of  some  modules.  If  a  complete  elimination  of 
the  feedbacks  is  not  possible,  they  are  reduced  in  number 
and  clustered.  An  example  of  a  fairly  large  Af -square  Matrix 
in  the  initial,  arbitrary  order  is  shown  in  Hg.  ISa  while  its 
clustered  state  is  shown  in  Fig.  ISb.  In  the  clustered  state 
the  system  is  hybrid — partially  hierarchic  and  partially  non- 
hieraichic.  A  software  tool  that  is  available  to  make  the  above 
transformation  is  described  in  ref.  9.  All  the  modules  in  one  of 
the  clusters  in  Fig.  tSb  may  be  regarded  as  a  new  supermodule, 
and  the  system  diagram  may  be  drawn  in  terms  of  these 
supermodules  as  shown  in  Fig.  16.  This  diagram  defines  a 
hierarchic  decomposition  of  a  system  because  the  data  flow 
from  the  top  of  the  pyramidal  hierarchy  to  the  bottom,  without 
reversing  the  flow  and  without  lateral  flow,  while  inside  of 
each  cluster  there  is  a  system  whose  modules  define  a  non- 
hietarchic  decomposition. 

The  Af-square  Matrix  structure  has  a  reflection  in  the  struc¬ 
ture  of  the  matrix  of  coefficients  in  eq.  S:  each  feedforward 
instance  in  the  former  gives  rise  to  a  Jacobian  matrix  located 
below  the  diagonal  in  the  latter  and  each  feedback  is  reflected 
in  a  Jacobian  above  the  diagonal.  Hence,  a  sequential  system 
without  feedbacks  has  a  matrix  of  coefficients  populated  only 
below  the  diagonal  so  that  eq.  5  may  be  solved  by  backsubsti- 
tution  of  the  right  hand  sides  without  factoring  of  the  matrix 
of  coefficients. 

SOP  Adapted  to  Hierarchic  System 
When  a  decomposed  system  has  a  hierarchic  structure,  its  SOP 
may  be  reorganized  to  include  separate  optimizations  in  each 
Mack  box.  This  SOP  version  was  introduced  in  ref.  10  and 
called  an  optimization  by  linear  decomposition.  It  has  found 
a  number  of  applications,  for  example,  it  was  the  basis  for  an 
algorithm  for  multilevel  structural  optimization  by  substructur- 
ing  in  ref.  II,  and  its  use  in  multidisciplinaiy  applications  was 
repotted  in  ref.  12  for  control-stnicture  interaction  and  in  ref. 
13  for  optimization  of  a  transpon  aircraft. 

Multilevel  optimization  of  a  hierarchic  system  by  a  linear  de¬ 
composition  exploits  the  top-down  Row  of  the  analysis  infor¬ 
mation.  At  the  bottom  level,  the  inputs  obtained  from  analysis 
at  the  iKxt  higher  level  and  the  appropriate  design  variaMes 
are  regarded  u  constanu  in  optimizadon  of  each,  bottom-level 
Mack  box.  Derivadves  of  each  such  opdmizadon  are  computed 
with  respect  to  these  input  consumts  by  means  of  an  algo¬ 


rithm  described  in  ref.  14  and  are  used  in  linear  extrapolations 
(hence  the  name  of  the  technique)  to  approximate  the  effect 
of  the  input  constants  on  the  optimization  results.  Optimiza¬ 
tions  in  the  black  boxes  at  the  next  higher  level  approximate 
their  influence  on  the  lower  level  optimizadon  by  means  of 
these  extrapolations.  Thus,  the  top  black  box  opdmizadon  is 
performed  taking  an  approximate  account  of  the  effect  of  its 
variables  (the  system  level  variables)  on  all  the  black  boxes  in 
the  hierarchic  pyramid.  As  mendoned  in  the  foregoing,  the  ad¬ 
vantages  of  the  SOP  exploiting  the  hierarchic  structure  of  the 
system  is  a  separadon  of  the  bottom  level  detailed  opdmiza- 
dons  from  the  ti^  level  system  optimization,  and  breaking  the 
large  system  optimization  problem  into  a  number  of  smaller 
opdmizadon  problems,  in  contrast  to  the  non-hierarchic  sys¬ 
tem  SOP  (Fig.  7)  in  which  optimizadon  is  performed  for  the 
system  as  a  whole.  However,  if  any  of  these  black  boxes 
in  a  hierarchic  system  contains  a  cluster  (see  discussion  of 
Fig.  16)  of  black  boxes  forming  a  non-hierarchic  system,  the 
non-hierachic  system  SOP  (Fig.  7)  may  be  used  to  optimize 
it  locally.  Hence,  both  methods  for  system  optimization  de¬ 
scribed  above,  the  one  based  on  the  linear  decomposition  (ref. 
10)  as  well  as  the  SOP  based  on  Fig.  7  flowchart  have  their 
place  in  optimization  of  a  general  case  of  a  hybrid  engineer¬ 
ing  system  that  exhibits  both  the  hierarchic  and  non-hierarchic 
structures  depicted  in  Fig.  16. 

As  reported  in  ref.  13,  the  linear  decomposition  method  was 
used  to  optimize  the  variables  of  configuration  geometry  and 
cross-sectional  structural  dimensions  of  a  transport  aircraft  il¬ 
lustrated  in  Fig.  17a  for  minimum  fuel  burned  in  a  prescribed 
mission,  under  constraints  drawn  from  the  disciplines  of  aero¬ 
dynamics,  performance  and  structures.  The  analysis  was  rel¬ 
atively  deep,  e.g.,  a  CFD  code  in  aerodynamics,  and  a  finite 
element  model  of  the  built-up  structure  of  the  airframe  struc¬ 
tures.  The  number  of  design  variables  was  over  1300,  and 
the  number  of  constraints  was  also  in  thousands.  Optimization 
was  conducted  decomposing  the  problem  into  a  three-level  hi¬ 
erarchic  system  shown  in  Fig.  17b.  A  sample  of  results  is 
depicted  in  Fig.  18  showing  a  smooth  convergence  of  the  fuel 
mass  and  the  structural  weight  in  only  4  to  6  cycles  (one  cycle 
comprised  the  top-down  analysis  and  the  bottom-up  optimiza¬ 
tions),  for  both  feasible  and  infeasible  initial  design. 

5.  GENERALIZATION  TO  ENTIRE  VEHICLE  DESIGN 
PROCESS 

The  approach  to  the  system  sensitivity  and  optimization  dis¬ 
cussed  in  the  foregoing  may  be  generalized  to  serve  the  entire 
design  process  as  shown  in  ref.  IS  using  as  an  example  a  def¬ 
inition  of  that  ^ocess  given  in  ref.  16.  The  process  defined 
in  ref.  16  is  a  conventional,  sequential  process  illustrated  in 
Hg.  19.  As  suggested  in  the  u|q)er  right  comer  of  the  flow¬ 
chart,  any  change  in  a  major  design  variable  such  as  the  wing 
or  engine  size  requires  reentry  into  the  sequence  and  repeti¬ 
tion  of  all  the  operations  in  the  chain.  However,  the  black 
boxes  foiming  the  sequence  are  also  forming  a  coupled  sys¬ 
tem  whose  diagram  is  depicted  in  Fig.  20.  The  arrows  in  the 
diagram  represent  the  data  flow  among  the  black  boxes,  exam¬ 
ples  of  the  data  being  defined  in  Table  2.  Application  of  the 
SSA  based  on  eq.  5  to  the  system  in  Fig.  20  leads  to  GSE  in 
the  format  shown  in  Fig.  21.  In  the  abbreviated  notation  used 
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in  that  figure,  Yij  stands  for  a  Jacobian  nutrix  Jij  defined  in 
eq.  3.  Solution  of  the  equations  shown  in  Fig.  21  yields  the 
SDD  values  that  answer  the  “what  if’  questions  implied  in  the 
upper  right  comer  of  the  flowchart  in  Fig.  19,  and  does  it  for 
all  the  variables  of  interest  simultaneously  and  without  repeat¬ 
ing  the  entire  chain  for  every  question.  The  SDD  values  may 
then  be  used  to  support  judgmental  design  decisions  and/or  to 
guide  a  formal  optimization  according  to  the  SOP  in  Fig.  7. 

fi.  CONCLUDING  REMARKS 

Design  of  an  engineering  system,  such  as  an  aircraft,  is  a 
formidable  task  involving  a  myriad  of  cross-influences  among 
the  engineering  disciplines  and  parts  of  the  system.  The 
time-honored  approach  to  that  task  is  to  decompose  it  into 
smaller,  more  manageable  tasks.  The  paper  outlines  some 
recently  developed  techniques  that  support  such  an  approach 
by  building  an  engineering  system  optimization  on  a  modular 
basis,  that  comprises  engineering  specialty  groups  and  their 
black  box  tools  and  allows  engineers  to  retain  responsibility  for 
their  domains  while  working  concurrently  on  marugeable  tasks 
and  communicating  with  each  other  by  means  of  sensitivity 
data.  The  modularity  and  concurrence  of  operations  map 
onto  the  familiar  structure  of  the  engineering  organizations 
and  are  compatible  with  the  emerging  computer  technology 
of  multiprocessor  computers  and  distributed  computing.  The 
only  major  new  requirement  is  the  generation  of  derivatives  of 
output  with  respect  to  input  in  each  specialty  domain. 

The  use  of  sensitivity  data  as  the  communication  medium  is  the 
distinguishing  feature  of  the  proposed  approach  and  represent 
a  major  improvement  over  the  present  practice  because  it  adds 
the  trend  information  to  the  function  value  information.  Both 
types  of  information  enhance  the  human  judgment  and  intuition 
while  being  readily  usable  in  guiding  the  formal  optimization 
procedures. 
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7.  Flowchart  of  the  System  Optimization  procedure  (SOP). 
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8.  Hypersonic  aircraft;  some  of  the  configuration  design 
variables. 
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10.  System  of  mathematical  models  for  hypersonic  interceptor  optimization. 
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17,  Optimization  of  a  transport  aircraft;  a)  Configuration:  b)  Hierarchic  system  of  modules. 
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18.  Sample  of  results  from  transport  aircraft  optimization. 
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19.  A  conventional,  sequential  design  process  for  aircraft. 
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20.  Black  boxes  from  Fig.  19  forming  a  system. 
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•  System  sensitivity  equations  of  design  represented  as  coupled  system 
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•  These  system  derivatives  answer  "What  if"  questions  regarding  these  variables 
without  reanalyzing  the  system 


21.  GSE  matrix  for  the  system  of  Fig.  20. 
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Abstract 

Formal  mathematical  optimization  methods  have  been 
developed  during  the  past  10  to  15  years  for  the  structural  de¬ 
sign  of  aircraft.  Together  with  reliable  analysis  programs  like 
finite  element  methods  they  provide  powerful  tools  for  the 
structural  design.  They  are  efficient  in  at  least  two  ways; 

•  producing  designs  that  meet  all  specified  requirements  at 
minimum  weight  in  one  step; 

•  relieving  the  engineer  from  a  time  consuming  search  for 
modifications  that  give  bener  results,  they  allow  more 
creative  design  modifications. 

MBB  has  developed  a  powerful  optimization  code  called 
MBB-Lagrange  which  uses  mathematical  programming  and 
gradients  to  fulfill  different  constraints  simultaneously  [  I ). 

Some  examples  depicting  the  successful  application  of 
the  MBB-LAGRANGE  code  are  presented.  Also  results  of 
other  optimization  codes  are  shown. 

The  paper  closes  with  an  outlook  how  the  optimization 
problem  could  be  enlarged  to  include  also  shape  and  size  of 
airplanes. 

IntmiltfaiOT 

To  improve  or  modify  a  design,  a  process,  a  procedure,  or 
any  given  task  into  a  "better"  direction,  is  refened  to  as  "opti¬ 
mization".  This  is  often  done  by  experience,  parametric  inve¬ 
stigations,  iterative  procedures,  by  experimental  testing  and 
modifications,  or  based  on  empirical  data.  Such  an  approach 
usually  leads  to  better  results  but  nobody  can  tell  how  far 
away  the  optimum  still  is  or  even  where  it  is.  A  more  efficient 
way  to  perform  this  task  is  provided  by  a  special  branch  of 
applied  mathematics,  called  optimization.  This  kind  of  opti¬ 
mization  changes  the  chosen  variables  in  a  design  problem  in 
a  way  to  achieve  the  best  value  for  an  objective  while  not  vio¬ 
lating  defined  constraints  that  represent  the  boundaries  of  the 
design  space. 

This  formal  optimization  was  rather  early  introduced  in 
economics  or  chemkal  engineering  due  to  the  linearity  of  the 
problems,  as  described  by  Ashley  in  an  excellent  overview 
paper  on  the  aeronautical  use  of  optimization  [2].  In  order  to 


use  the  potential  of  mathematical  optimization,  it  is  necessary 
to  describe  the  physical  nature  of  the  problem  in  a  way  that  all¬ 
ows  the  use  of  optimization  algorithms. 

In  structural  design,  finite  element  methods  together  with 
modem  computers  have  provided  tools  that  allow  to  analyse 
complex  structures  with  high  accuracy.  These  were  main  essen¬ 
tials  to  initiate  development  and  application  of  optimization 
programs  for  structural  design  in  1970.  Approximately  at  the 
same  time,  composite  materials  were  introduced  in  aerospace 
design.  They  offer  an  infinite  variety  to  combine  their  highly 
anisotropic  elastic  properties  for  any  specific  combination  of 
design  requirements.  For  a  more  efficient  use  of  these  mate¬ 
rials,  optimization  programs  are  required  to  handle  the  comple¬ 
xity  of  the  problem,  especially  if  additional  requirements  besi¬ 
des  strength  are  involved  in  the  problem  [3J.  During  the  last 
decade  considerable  effwt  has  been  spent  to  develop  modem 
structural  optimization  procedures,  using  efficient  mathematical 
optimization  algorithms  as  well  as  optimality  criteria  which  sa¬ 
tisfy  all  requirements  simultanously  and  find  optimal  values  of 
the  design  variables  by  direct  computation.  The  increasing  em¬ 
phasis  of  aeroelastic  considerations  is  shown  in  Fig.  1,  which 
was  taken  from  [4]. 
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Fig.  1  EVOLUTION  OF  AEROELASTIC  CON- 
SIDERATIONS  IN  FIGHTER  AIRCRAFT 
DESIGN 


STRUCTURAL  OPTIMIZATION  IN  THE  GENERAL 

DATAF 

LOW 


The  use  of  structural  optimization  tools  during  the  preli¬ 
minary  design  stage  of  an  advanced  aircraft  gives  the  follo¬ 
wing  potential  improvements: 

•  satisfies  the  requirements  of  new  aircrafts 

•  minimizes  the  objective  (weight) 

•  increases  the  quality  of  products 

•  shortens  the  development  phase 

•  increases  chances  of  the  company  in  competition. 


In  order  to  do  this,  an  appropriate  mathematical  pro¬ 
gramming  procedure  has  to  be  embedded  in  the  general  data 
flow,  which  is  depicted  in  Fig.  2  and  Fig.  3  taken  from  (5). 

For  the  general  data  flow  I  refer  to  ( 1 1. 


Fig.  2  EXTERNAL  GEOMETRY  IN  DATA 
FLOW 


Fig.  3  NEW  TECHNOLOGIES  OF  RECENT 
AIRCRAFT 


STRUCTURAL  OPTIMIZATION  AT  MBB 

MBB  has  developed  its  own  structural  optimization  sy¬ 
stem  called 

MBB-LAGRANGE 

The  performance  and  requirements/constraims  of  this 
new  program  system  are  shown  below: 

•  Requirements 

Finite  Element  Structure 

•  Structure  Variables 

Skin  Thickness 
Balance  Masses 
Fibre  Directions 
Grid  Point  Coordinates 


These  figures  show  a  typical  flow  of  geometric,  aerody¬ 
namic,  structural  and  other  data  which  arc  used  during  the  de¬ 
sign  phase  of  an  aircraft.  The  improved  productivity  is  a  re¬ 
sult  of  the  integrating  effects  of  the  structural  optimization. 
Shorter  time  of  development  is  realised  and  fewer  data  trans¬ 
fers  go  wrong. 


At  the  present  time  the  development  of  new  airplanes  is 
influenced  by  new  techniques,  such  as  flutter  suppression,, 
CCV-configuration,  gust  load  alleviation  etc.  (Fig.  3).  In  ad¬ 
dition  to  stress,  displacennent,  aeroelastic  and  dynamic  con¬ 
straints  an  integrated  design  involves  all  these  techniques  and 
the  optimization  procedures  must  be  extended  for  these  new 
constraints.  A  reliable  cqjtimization  code  is  the  basis,  which 
allows  parametric  investigations  and  weight  penalties  to  be 
evaluated  properly. 


Constnunts 

MinVMax.  -  Variable 

Stresses 

Strains 

-  Deformations 
Flutter  Speed 
Divergence  Speed 
Aeroelastic  Fificiencies 
Eigen  Frequencies 
Element  Stability 
Dynamic  Response 
Weight 

For  mote  information  reference  is  made  to  [1]. 
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•  Muitiobiective  Function  r(x)  =»  Min. 

Vector  Optimization  =  'Trade  Off'  Studies  of  Convex 

Combination  of  Objectives 

The  program  architecture  is  organized  according  to  the 
concept  of  H.  A.  Eschenauer  [6]  with  the  main  parts  optimi¬ 
zation  algotitiim,  optimization  model  and  strtictural  analysis 
including  sensitivity  analysis. 

The  corresponding  optimization  models  are  based  on  the 
general  nonlinear  programming  problem  according  to  [6], 
The  design  variables  x  are  cross  sectional  areas  of  trusses  and 
beams,  wall  thicknesses  of  membrane  and  shell  elements,  la¬ 
minate  thicknesses  for  every  single  layer  in  composite  ele¬ 
ments  or  nodal  coordinates  for  geometry  optimization  pro¬ 
blems.  The  constraints  in  form  of  inequalities  may  be  any 
combination  of  displacements,  stresses,  strains,  buckling 
loads,  aeroelastic  efficiencies,  flutter  speed,  divergence 
speed,  natural  frequencies,  dynamic  response  and  design  va¬ 
riables  (7]. 

In  the  case  of  scalar  optimization,  the  objective  function 
f(x)  often  includes  the  structural  weight  or  another  linear 
combination  of  the  design  variables.  However,  it  is  also  pos¬ 
sible  to  defme  one  of  the  constraint  functions  as  objective  and 
to  introduce  the  weight  as  constraint  at  the  same  time.  If  vec¬ 
tor  optimization  problems  are  under  consideration,  then  opti¬ 
mization  strategies  p[f(s)l  according  to  [6]  ensure  the  trans¬ 
formation  to  scalar  substitute  problems. 

It  is  necessary  to  provide  several  different 
algorithms,  because  there  is  no  known  single  algorithm  which 
is  adapted  to  every  type  of  problem.  Some  of  the  algorithms 
which  are  implemented  in  LAGRANGE  are  shown  below: 


•  IBF 

Inverse  Barrier  Function, 

•  MOM 

Method  of  Multipliers, 

•  SLP 

Sequential  Linear  Programming, 

•  SRM 

Stress  Ratio  Method, 

•  R0P1,RQP2 

Recursive  (Quadratic  Programming 

•  GRG 

Generalized  Reduced  Gradients 

A  more  extensive  explanation  is  given  in  (1|. 

Fig.  S  denoonstrates  the  iteration  steps  to  achieve  the  desi¬ 
red  flutter  speed  of  900  lets  starting  at  700  kis  for  the  initial 
fully  stressed  design.  Results  of  flutter  calculations  for  botii 
designs  are  plotted  in  Fig.  6.  It  should  be  nKntioned  also  that 


Fig.  4  FLUTTER  SENSITTVITIES  FOR 
WING  COVER  SKIN 


Fig.  5  RESULTS  OF  REDESIGN  STUDY 


Early  Investigations  with  FASTOP 

In  1979  the  structural  optimization  program  FASTOP 
(Flutter  and  Strength  Optimization  Program)  was  acquired  by 
MBB.  The  capability  of  the  program  was  extended  to  be  able 
to  analyse  also  the  static  aeroelastic  behaviour  of  structures. 
Main  features  of  the  program  were 

•  a  fully  stressed  design  for  static  requirements  that  usually 
goes  near  optimal  results. 

•  a  flutter  redesign  for  a  defined  minimum  speed. 

The  program  was  extensively  used  for  wing  design  stu¬ 
dies  (8).  An  example  for  the  flutter  redesign  of  a  stress  desig¬ 
ned  (optimised)  wing  is  shown  in  Fig.  4.  Essential  elements 
that  are  sensitive  for  flutter  are  marked  for  the  upper  and  lo¬ 
wer  cover  skin.  At  that  time  the  program  could  not  handle  the 
different  layers  of  composite  materials  individually  as  design 
variables.  The  design  space  was  limited  by  predetermined  fi¬ 
bre  orientations  and  the  proportion  of  the  individual  layer 
thickness  in  the  total  laminate. 


Fig.  6  FLUTTER  ANALYSES 
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AEROELASTIC  PROBLEMS  AND  STRUCTURAL 
DESIGN  OF  A  TAILLESS  CFC-SAILPLNE  f91 


Tailless  planes,  sometiines  also  called  "Flying  Wings", 
have  always  been  a  challenge  in  airplane  design.  They  offer  a 
great  potential  in  perfonnance  compared  to  conventional  de¬ 
signs  because  of  less  surface  (parasite  drag),  less  weight,  and 
less  trim  drag.  Although  this  has  been  known  for  a  long  time, 
tailless  planes  never  have  experienced  the  success  one  might 
expect.  Many  carefully  designed  tailless  gliders  had  to  be  re¬ 
designed  or  designed  completely  new  after  first  flight  tests 
because  of  strange  instabilities,  which  were  often  not  under¬ 
stood  or  misinterpreted. 

It  was  an  interesting  task  for  the  "Akademische  Flieger- 
gruppe  Braunschweig"  to  launch  a  tailless  ailplane  project  for 
the  IS  meters  standard  class  in  1983.  During  flight  test  with  a 
1/3  scaled  remotely  piloted  model  a  severe  instability  occured 
at  very  low  speed.  Flutter  calculations  using  data  from  a 
ground  resonance  test  showed  that  coupling  of  the  rigid  body 
shon  period  mode  with  the  first  elastic  mode  caused  the  phe¬ 
nomenon.  Solving  this  problem  is  a  multidisciplinary  task. 
MBB  offered  assistance  to  redesign  the  wing  with  the  help  of 
modem  optimization  programs. 

By  applying  these  codes  the  flutter  speed  could  be  increa¬ 
sed  to  an  acceptable  level  with  small  modifications  of  the 
wing  root  geoittetry,  a  new  design  of  the  main  spar,  and  by 
the  use  of  a  new  high  modulus  fiber  type.  With  a  small  weight 
penalty-compared  to  the  initial  design  -  the  flutter  speed  could 
be  doubled. 

Sailplanes  have  achieved  a  very  high  technological  level 
during  the  last  20  years,  mainly  due  to  fiber  composite  struc¬ 
tures  and  improved  aerodynamics. 

Further  improvements  can  be  expected  only  from  small 
detail  modification  or  expensive  projects  like  variable  wing 
geometry.  For  this  reason  an  unconventional  design  concept 
like  the  tailless  wing  is  a  challenge  for  designers.  It  offers  se¬ 
veral  advantages  like 

•  less  parasite  drag 

•  less  weight 

•  less  construction  effort 

due  to  the  missing  rear  fuselage  and  the  tail.  If  the  airfoil  is 
designed  carefully  for  a  small  pitching  moment,  the  flying 
wing  will  not  have  a  higher  profile  drag  compared  to  conven¬ 
tional  planes.  If  the  vertical  tail  is  integrated  in  winglets,  the 
advantage  of  less  induced  drag  can  be  explored  without  addi¬ 
tional  weight  penalties. 

The  SB  13  project.  Fig.  7,  shows  performance  improve¬ 
ments  of  up  to  10%  conqiaied  to  existing  conopetitors  in  the 
15-meters  standard  class,  as  indicated  in  Fig.  8  for  the  veloci¬ 
ty  polar.  Table  1  gives  some  main  design  parameters. 


Fig.  7  3-SIDE  VIEW  OF  THE  TAILLESS 
SAILPLANE  PROJECT  SB13 

HorizontaY  Spe«d  kia/h 


SBI3 


Wine 

Span 

15  m 

Area 

11.6  m2 

Aspect  Ratio 

19.4 

Dihedral 

+  4” 

Tvrist 

- 1.5°  outboard 

Wing  Section 

HQ  34  N/14.83  inboard 

HQ  36  N/15.12  outboard 

Winrdets 

Length 

1.25  m 

Area 

0.675  m^ 

Aspect  Ratio 

2.31 

Pn^e 

FX-71-L  15(V30 

ESKtyS 

Length 

3.02  m 

Width 

0.66  m 

Height 

0.84  m 

Landing  gear 

2  retractable  wheels,  spring- 
suspended 
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Weights 

Weight  empty 

2240  N 

Payload 

700-  llOON 

Water  ballast 

max.  1330  N 

Gross  weight 

2940  -4270  N  , 

Wing  Loading 

248  -  360N/m" 

Performance 

V  ■ 
mm 
mm 

70kni/h 

max 

mm.  sink  speed 

210km/h 

0.53  m/s 

max.  L/D-ratio 

43.5 : 1 

TABLE  1 

TECHNICAL  DATA 
PLANE  SB13 

Although  tailless  aircraft  have  been  studied  almost  since 
the  beginning  of  aviation  they  have  never  experienced  the 
success  one  might  expect.  One  reason  for  the  lack  of  success 
is  described  in  [10]  and  [11]  as  the  extreme  difficulty  of 
achieving  satisfactory  unaugmented  handling  qualities,  con¬ 
trol  and  dynamic  stability. 

Probably  the  most  experienced  designer  of  tailless  planes 
was  A.  Lippisch  with  numerous  powered  and  unpowered  de¬ 
signs  [12].  He  reported  about  several  difficulties  and  crashes, 
caused  among  others  by  longitudinal  oscillations  or  "unsatis¬ 
factory  handling  qualities"-  The  Honen  brothers  also  desig¬ 
ned,  constructed,  and  tested  various  tailless  planes  between 
1936  and  1960  [13]. 

Among  the  few  successful!  tailless  sailplane  were  the 
single-seat  AV-36  and  the  twin-seat  AV-22  by  Charles  Fau- 
vel  [14],  Fig.  9,  and  the  very  similar  looking  designs  from  J. 
Marske  [15]. 


Fig.  9  FAUVEL  TAILLESS  SAILPLANES 
AV.22  AND  AV.36 


To  snidy  stability  and  handling  qualities  of  the  SB  13,  a 
remotely  piloted  1/3  scale  model  was  built  and  flown.  The 
handling  qualities  showed  no  problems,  but  an  unexpected  in¬ 
stability  in  the  longitudinal  motion  occured  at  very  low 
speeds.  When  a  ground  resonance  test  and  a  flutter  calcula¬ 
tion  was  performed,  it  could  be  shown  that  the  reason  for  the 
instability  was  the  coupling  between  the  rigid  body  short  pe¬ 
riod  mode  and  the  first  symmetric  structural  mode. 

After  the  problem  was  solved  analytically  by  the  means 
of  aeroelastic  tailoring  and  the  a|q>Iication  of  a  new  carbon  fi¬ 
ber,  a  paper  about  a  very  similar  design  study  at  Cranfield 


was  published  [16].  This  project,  called  "Ricochet"  showed 
the  same  aeroelastic  behaviour  as  the  SB  13.  Fig.  10  shows 
the  design  and  gives  some  design  parameters.  Because  the 
flutter  problem  could  not  be  solved  in  this  case,  the  project 
was  finally  given  up.  But  this  study  is  the  first  one  known  to 
the  authCHS  which  identified  the  problem  craiectly. 


Paraaeter 

The  Ricochet 

Hftterial  of  construction 

AluminiuD  alloy 
(6061-T6) 

Span 

15m 

Wing  area 

10.26 

Aspect  ratio 

22.93 

Wing  root  chord 

0.73  B 

Ving  tip  chord 

0.50  ■ 

Sweep  angle 

13° 

Mass  of  each  wing 

50  Kg 

Puselage  mass  with  equipments 

65  Kg 

Fig.  10  RICOCHET  SAILPLANE  PROJECT 

Flutter  Calculations  for  the  RPV-Modei 

To  investigate  the  flutter  behaviour  more  thoroughly,  a 
ground  resonance  test  was  petfotmed  at  the  DFVLR,  Institute 
for  Aeroelasticity  in  Gbttingen.  Fig.  1 1  shows  the  test  instal¬ 
lation,  Table  2  gives  the  main  results  for  two  configurations, 
where  configuration  II  contains  additional  fuselage  mass  for 
non-structural  items.  The  first  bending  mode  for  both  confi¬ 
gurations  is  shown  in  Fig.  12.  Several  flutter  calculations 
were  performed  using  the  described  data.  If  the  rigid  body 
tnodes  are  ignored,  the  first  structural  mode  shows  divergence 
in  the  flutter  calculation  as  indicated  in  Fig.  13  for  free-free 
boundary  condition. 


Fig.  11  GROUND  RESONANCE  TEST 


EQUIPMENT  FOR  SCALED  RPV 
MODEL 


FREQUENCY  [fiz] 
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tot 


By  (k9M*I  0.683 
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Fig.  14  gives  a  isometric  view  of  mode  1  and  Fig.  IS 
shows  time  clearly  the  local  deformations  at  the  quarter 
chord  line  for  both  configuration. 

In  Fig.  16  the  eigenvalues  of  the  shcnt  period  are  plotted 
vs.  airspeed  for  the  full  scale  SB13  and  die  conventional 
SBl  1.  The  frequency  of  the  SB13  is  almost  three  times  that  of 
the  SB  11  while  the  danqiing  is  smaller.  This  difference  is 
mainly  due  to  the  small  pitch  moment  of  inertia  of  the  SB13. 
Table  3  gives  a  comparison  of  inqxmant  parameters  for  the 
longitudinal  motion  of  the  model,  the  SB13  and  the  SBl  1. 

The  equations  for  the  short  period  mode  are 


for  the  frequency,  and 


for  the  damping.  Here  C^a  and  are  functio^rf  geome¬ 
try  only,  Cw^  also  depends  on  the  c.g.  location.**f(fc|,musi 
be  identicalfor  the  model  and  the  full  scale  version,  if  geome 
trical  proportions  are  similar  and  the  static  longitudinal  stabi¬ 
lity  is  equivalent  The  last  parameter  for  comparison  is  the  re¬ 
lative  mass  density. 


>'1 


2m 


(3) 


As  table  3  shows,  these  terms  are  identical.  The  construction 
of  the  model  is  similar  to  the  one  of  a  modem  sailplane.  The¬ 
refore  it  can  be  expected  that  there  is  also  elastic  similarity 
(replica). 

With  the  results  from  the  ground  resonance  test  the  rigid 
body  coupling  with  the  elastic  mode  results  in  a  flutter  speed 
of  53  kmAi  for  configuration  I,  Fig.  17,  and  due  to  the  smaller 
fj  only  44  km/h  for  configuration  11,  Fig.  18.  If  we  assume  li¬ 
nearity  between  flutter  speed  and  first  elastic  ttrade  frequency, 
we  get 


fp 


(4) 


(L  =  large  scale,  M  0  model)  for  the  velocity  scale.  With  the 
length  scale  X =  2.71  we  got 

.  «.3  “ 

for  the  SB13.  If  we  demand  a  flutter  speed  equal  to  the  maxi¬ 
mum  velocity  Vjj  «  283  km/h  for  the  SB13,  this  would  requi¬ 
re  a  first  structural  mode  frequency  of  6.7  Hz. 


SB  13-Model 

9  1] 

SB  1 1  •  IS  « 

Confttiratto*  11 

Fontcrflaet  1* 

»  fkal 

<Z.( 

200 

465 

o,as4 

113 

835 

1,  1«1 

0.26S 

Q.B3S 

1,340 

s  !••) 

1.33 

11,0 

13.32 

c  Ul 

0.271 

0.73S 

0.772 

“i 

S7.07 

M.M 

73.83 

S.72 

S.72 

S.7I 

'«<t 

-Z.M* 

-2.749 

-20.49 

0.194 

0,194 

0,649 

iCp 

-0.10 

-O.IO 

•0.10 

Table  3  LONGITUDINAL  MOTION  PARA¬ 
METERS  FOR  SB13,  AND  A  COM¬ 
PARABLE  EXISTING  SAILPLANE 
(SBll) 


Fig.  17  FLUTTER  CALCULATION  WITH 
RIGID  BODY  DEGREES  OF  FREE¬ 
DOM  FOR  MODEL  CONnGURA- 
TIONI 


I 

I 


Frequency  Damping  L»J 
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Fig.  18  FLUl  TER  CALCULATION  FOR 
CONFIGURATION  II  INCLUDING 
ADDl  nONAL  MASS  FOR  THE 
PILOT 


POBsible  Solutioite  to  Increase  Uie  Flutter  Speed 

As  described  F  fore,  the  flutter  problem  is  caused  by  the 
coupling  of  strucpiml  bending  mode  (Bl)  and  rigid  body 
shon  period  mode  iSl).  Obviously  a  separation  of  uie  two 
frequencies  would  be  favourable. 

•  To  reduce  mode  SI  frequency  one  must  increase  the  pitch 
moment  of  inertia  largely  which  is  not  possible  with  a 
tailless  configuration.  As  also  the  Ricochet  study  showed, 
this  is  the  only  imponant  parameter  for  the  shon  period 
mode.  Static  longitudinal  stability,  or  wing  sweep  angle 
do  not  improve  the  rigid  body  motion.  Fig.  19  shows  the 
change  in  flutter  speed  with  the  pitch  moment  of  inertia 
for  the  Ricochet  s^  plane. 

•  Changing  the  configuration  in  a  way  that  the  first  elastic 
mode  shape  (SI)  frequency  shows  no  reduction  with  in¬ 
creasing  airspeed  would  require  a  congrleiely  new  design. 
Because  a  large  effort  had  already  been  invested  in  the 
aerodynamic  layout,  this  solution  was  not  desirable. 

•  Active  Control  Technology  is  a  good  method  to  extend 
the  flight  envelope.  Wykes  [17]  proposed  a  CSAS  sy¬ 


stem  for  the  U.S.  forvyard  swqrt  wing  fighter  project,  where 
the  same  coupling  between  rigid  body  and  elastic  structure 
occurs.  Although  tested  for  several  military  and  commer¬ 
cial  aircraft  projects  sucessfully,  ACT  is  not  a  feasible  solu¬ 
tion  for  sailplanes.  It  would  require  power  supply  and  a 
complicaied  sensor-,  control-,  and  actuation  system. 

•  To  change  the  aeroelastic  behaviour  using  mass  balance  by 
addition  of  lumped  masses  does  not  improve  the  situation 
with  feasible  arrangements. 

•  Decoupling  the  pilot  ftom  the  fuselage  to  change  the  criti¬ 
cal  mode  sfaape^equency  would  result  in  an  unfavourable 
sensing  system  for  the  pilot  In  addition  all  spring  systems 
for  this  purpose  would  have  large  anqrlitudes  (non  lineari¬ 
ties)  under  load,  or  would  require  too  much  volume  (air 
bag). 

•  The  only  practicable  solution  is  a  combination  of  a  strucm- 
ral  redesign  (with  small  modifications  in  the  wing  root  geo¬ 
metry),  using  high  elastic  modulus  carbon  fibers  to  increase 
the  first  elastic  frequency  and  tailor  the  wing  fw  a  different 
aeroelastic  behaviour  (exploiting  the  anisotropic  material 
properties  to  change  mode  shapes).  This  procedure  -  finally 
selected  -  will  be  decribed  in  the  next  chapter. 


Fig.  19  FLUTTER  SPEED  VS,  PITCH  MO¬ 
MENT  OF  INERTIA  FOR  THE 
RICOCHET  PROJECT 


The  first  handic^  in  the  application  of  TSO  and  FASTOP 
for  "aeroelastic  tailoring"  the  wing  was  caused  by  the  lack  of  ri¬ 
gid  body  modes  in  both  computer  programs.  Therefore  we  had 
to  choose  substitution  systems  to  describe  the  critical  flight  me¬ 
chanical  mode.  This  can  be  achieved  by  defining  soft  springs 
between  the  structure  and  an  earthed  poinL  The  softness  of  the¬ 
se  springs  is  limited  by  numerical  proMems  in  the  stiffness  ma- 
ffix.  Unfortunately  this  system  caused  other  problems.  The  very 
flexible  wing  attachment  did  not  allow  to  use  the  great  advanta¬ 
ge  of  TSO,  the  simultaneous  optimization  for  different  objecti¬ 
ves  and  constraints.  The  soft  attachment  caused  too  high  deflec¬ 
tions  under  load  for  the  strength  design. 

Only  limited  potential  of  aeroelastic  tailoring  sailplane- 
wings  is  available,  consmined  by  the  extremely  high  aspect  ra¬ 
tio  and  slenderness  of  these  wings.  The  main  fiber  direction  can 
be  swept  only  within  small  litttits.  Additionally  there  is  only  a 
small  number  of  45  °  ■  plies  which  makes  it  itttpossible  to  use 
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them  to  change  the  elastic  behaviour. 

Usually,  modetn  sailplane  wings  are  designed  as  shown 
in  Hg.  20.  There  is  one  main  spar  with  the  fluge  fabricated 
from  unidirectional  rovings.  The  unsional  forces  are  carried 
by  ±  45°  plies  in  a  sandvrich  shell  constructioiL  Alternatively 
a  shell  construction  with  coupled  bending  and  torsion  plies 
was  investigated  first.  Due  to  the  small  number  of  required 
plies  with  a  still  very  small  box  chord,  this  design  was  given 
up  later,  because  it  showed  no  improvements  and  is  also  very 
difficult  to  make  manually. 

The  final  solution  has  a  two  web  main  spar  (with  0°-phes) 
and  an  uncoupled  torsion  shell  (±  45°  -plies),  described  in  de¬ 
tail  later  on. 


Fig.  20  TYPICAL  WING  SECTION  FOR 
MODERN  SAILPLANES 

TSO-Calculations: 

This  program  describes  the  wing  structure  as  a  plate  mo¬ 
del.  Therefore,  check  calculations  were  necessary  to  see  if  the 
plate  theory  can  be  used  for  slender  wings.  For  a  constant 
chord,  constant  thickness  beam  with  an  afreet  ratio  of  20,  the 
TSO  results  could  be  confirmed  with  analytical  beam  themy. 
Hg.  21  shows  the  idealization  of  the  wing  box  within  the 
planform  geometry.  This  plot  also  shows  typical  thickness 
contours  for  the  bending  layers. 


Fig.  21  WING  IDEALISATION  FOR  TSO 
PROGRAM 


As  mentioned  befoie,  the  posribilities  for  aeroelastic  tailo¬ 
ring  were  limited.  Due  to  gemnetrical  constraints,  a  main  fiber 
sweep  angle  of  3°  forward  was  the  maximum.  Because  of  the¬ 
se  limitations,  several  different  materials  and  material  combi¬ 
nations  were  used  very  early.  Table  4  gives  some  typical  re¬ 
sults  from  these  optimization  runs. 


i  Carbon  Fibers  I 

Glass 

Fiber 

High 

Tension 

High 

Modulus 

Fiber  Volune  Ratio 

0.4 

0.65 

0.52 

E,,  Ifi/mf] 

30348 

131835 

195103 

^22 

6481 

7669 

6089 

G 

2120 

3220 

3434 

0.310 

0.3115 

0.272 

s  /s/emy 

1.696 

1.510 

1.469 

0.0129 

0.00455 

0.001287 

^22 

0.00347 

0.0022 

0.00266 

0.0188< 

0.01242 

0.01165 

Table  5  MATERIAL  PROPERTIES  USED 
IN  THIS  STUDY 


Run-No. 

Description 

Main  Fiber 
sweep  Angle 
(positive  fWd.) 

0'-p1ie$ 

/ksi/lj 

m 

initial  Design 

0° 

HT 

HT 

2.42 

74.1 

60.2 

increased  bending 
stiffness  with  HH-Fibers 

0° 

HM 

HT 

3.14 

99.0 

61.1 

additional  0*-p1ies 

0° 

HM 

HT 

4.31 

140.8 

79.7 

swept  0*-p1ies  skin  thick¬ 
ness  as  no.  2 

♦  2.5* 

■ 

HT 

3.40 

111.5 

61.1 

swept,  optiaized  lay  up 

♦  2.5* 

HT 

3.76 

115.0 

64.5 

*  3.0°,  linear  thickness 
distribution 

♦  3.0* 

U 

HT 

4.76 

162.2 

79.6 

B 

free  thickness  distribution 
lest  weight,  higher  than  6 

♦  3.0° 

HT 

4.8S 

178.5 

79.2 

8 

*  S.O,  glass  fiber  torsion 
plies,  sweep  angle  not  practl- 
calbe 

♦  5.0* 

B 

glass 

fibers 

187.6 

79.8 

Table  4  TSO  CALCULATIONS  FOR  SB13 
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RUN  1  (Initial  Desion) 
ROOT 


f,  =  2.42  Hz 


Fig.  22  TSO  MODE  SHAPE  1  FOR  INITIAL 
AND  OPTIMIZED  DESIGN 


Rather  soon,  it  became  obvious  that  the  flutter  problem 
could  not  be  solved  with  presently  used  materials. 

The  advent  of  high  modulus  carbon  fibers  provides  a 
Young's  modulus  30%  higher  than  in  presendy  used  high  ten¬ 
sion  fiber.  Table  3  gives  a  comparison  of  material  properties 
for  unidirecdonal  layers.  Fig.  22  shows  the  first  elastic  mode 
shape  for  the  inidal  design  and  for  the  high  modulus  fiber 
with  swept  0°-direcnon. 

In  the  flutter  calculation  results  (Fig.  23)  for  the  optimi¬ 
zed  design  No.  8  the  improved  flutter  behaviour  is  shown.  It 
should  also  be  mentioned  that  the  rigid  body  damping  (mode 
1)  is  predicted  wrongly. 


Fig.  23  TSO  FLUTTER  CALCULATION  RE¬ 
SULTS  FOR  OPTIMIZED  DESIGN 

FASTOP-CalgiilatInnn; 


Fig.  24  shows  the  finite  element  idealization  for  the  FA- 


STOP  model.  The  skin  is  described  by  membrane  elements 
whereas  the  ribs  and  spars  are  modelled  with  shear  panels.  Al¬ 
though  FASTOP  can  consider  fiee-fiee  conditions  in  the  vi¬ 
bration  analysis,  it  is  not  capable,  to  handle  rigid  body  modes 
separately.  They  can  only  be  superposed  to  the  elastic  modes. 
For  this  reason,  it  was  necessary  to  simulate  this  mode  with 
s(^  spring  attachments  and  cantilever  conditions  in  the  vibra¬ 
tion  analysis. 


Fig.  24  FIRST  STRUCTURAL  MODEL  FOR 
FASTOP  CALCULATIONS 


Instead  of  sweeping  the  fibers  within  the  spar  and  manu¬ 
facturing  the  spar  with  prepregs,  a  new  model  with  a  swept 
spar  inside  the  wing,  fabricated  conventionally  bom  rovings, 
showed  better  results  in  the  flutter  behaviour.  Fig.  23  shows 
the  new  idealization.  To  allow  a  higher  sweep  angle  for  the 
spar,  the  wing  planform  was  modified  in  the  inboard  section. 
The  leading  edge  sweep  angle  is  reduced  with  three  kinks,  the 
trailing  edge  sweep  angle  is  reduced  with  two  kinks.  The  in¬ 
tention  was  to  bring  the  main  spar  closer  to  the  pilot’s  mass 
without  increasing  the  wing  area. 


Fig.  25  MODIFIED  STRUCTURAL  MODEL 


Old  Design 


Redesign 


Fig.  26  WING  GEOMETRY  MODinCA- 

TIONS  FOR  IMPROVED  FLUTTER 
BEHAVIOUR 
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Fig.  26  shows  these  modifications. 

Using  HM-instead  of  HT-fibers  for  the  unswept  spar  in¬ 
creased  the  flutter  speed  from  1 10  km/h  to  210  km/h  (+  90%), 
sweeping  the  spar  3  degrees  gives  122  km/h  for  HT-fibers  (+ 
11%),  and  237  km/h  for  HM-fibers.  If  HT-caibon  fibers  are 
used  for  the  torsion  layers  instead  of  glass  fibers,  the  flutter 
speed  is  3.5  %  higher.  Table  6  gives  a  summary  of  these  re¬ 
sults. 

Fig.  27  shows  the  mode  shapes  and  flutter  calculation  re¬ 
sults  for  the  best  design.  This  wing  has  a  weight  of  67.7  kg 
compared  to  60.0  kg  for  the  initial  design  (-i-  12.8  %)  but  the 
flutter  speed  is  1 14  %  higher  !  Further  calculations  were  ne¬ 
cessary  for  different  water  tank  positions  in  the  wing.  Fig.  28 
shows  two  possible  solutions  which  do  not  decrease  the  flut¬ 
ter  speed.  Because  the  water  is  positioned  very  close  to  the 
nodal  line  of  mode  1,  the  first  ftequency  does  not  drop  more 
than  6%  while  the  shon  period  mode  is  more  than  10%  smal¬ 
ler  due  to  the  higher  moment  of  inertia.  The  flutter  calculation 
for  water  ballast  configuration  II  (No.  7  in  Table  6)  is  shown 
in  Fig.  29. 

Flutter  calculations  for  antisymmetric  modes  were  also 
performed.  Because  the  first  mode  is  higher  than  6  Hz  in  this 
case,  there  is  no  coupling  with  low  frequency  modes.  Higher 
modes  are  separated  without  tendency  to  couple  up  to  400 
km/h. 

The  final  configuration  was  analysed  using  a  different  ap¬ 
proach  whereby  the  rigid  body  mode  frequencies  could  be 
described  more  accurately  (0.01  Hz  for  the  z-translation  and 
0.05  Hz  for  the  rigid  rotation  at  0  airspeed).  This  influence 
improved  the  flutter  speed  considerably.  Table  6  gives  also 
the  results  using  this  method. 


From  these  results  it  could  be  expected  that  the  flutter 
speed  will  be  sufficiently  high  to  clear  the  full  flight  envelope 
up  to  the  maximum  speed  =  210  km^  (including  a  safety 
margin).  The  predictions  were  verified  during  flight  test  which 
happened  in  the  year  of  1985.  The  airplane  has  been  flown 
ever  since  and  has  never  shown  any  structural  instability. 

In  summary  it  can  be  said  that  most  tailless  sailplanes 
seem  to  have  the  great  disadvantage  of  high  frequency  shon 
period  modes  compared  to  conventional  constructions.  To  pre¬ 
vent  flutter  because  of  coupling  with  the  first  elastic  mode,  the 
wing  has  to  be  fabricated  from  extremely  stiff  materials.  New 
fibers  with  a  very  high  elastic  modulus  could  provide  the  ne¬ 
cessary  stiffiiess  for  the  SB  13  wing. 

To  overcome  these  difficulties  easier  in  other  tailless  airc¬ 
raft  designs,  it  could  be  possible  that  sweeping  the  wing  for¬ 
ward  and  thus  having  the  pilot  in  front  of  the  wing,  might 
change  the  first  elastic  mode  in  a  way  that  the  coupling  with 
the  shcMt  period  mode  will  be  delayed  to  higher  speeds.  J.  J. 
Marske  gathered  a  lot  of  experience  in  die  design  and  con¬ 
struction  of  several  tailless  sailplane  [15].  His  final  solution 
was  a  swept  forward  configuration  which  had  no  stability  pro¬ 
blems  and  showed  good  handling  and  perframance  characteri¬ 
stics.  For  the  SB  13  design,  a  swept  forward  solution  was  not 
possible  because  of  the  slender  fuselage  with  a  carry-through 
main  spar.  And  it  is  not  possible  if  the  winglets  are  used  as 
vertical  tails. 


Run- 

(Ho,  Description 

Main  Fiber 
Direction 
(pos.  fwd.) 

Ma 

O'-plies 

terials 

♦  45°-pl1es 

"wing 

CkJ 

per  side 

free/ fret 
divergence 
speed 
Ocm/hJ 

Flutter 

FASTOP 

Speed  Tltm/lj 
Check -pgm. 

1 

Initial  De¬ 
sign 

0” 

HT 

HT 

2.59 

2.74 

60.0 

186.0 

110.6 

120.0 

2 

HM-fibers 

0° 

HM 

HT 

5.27 

5.19 

69.5 

342.0 

210.2 

236.0 

3 

modified  main 
spar  HT-fibers 

3° 

HT 

HT 

2.85 

2.93 

60.0 

207.0 

122.4 

140.0 

4 

modified  spar 
HH  *  glass 
fibers 

3° 

HM 

glass 

fibers 

5.64 

5.46 

73.4 

390.0 

229.0 

260.0 

5 

HM  ♦  HT 
fibers 

3' 

HM 

HT 

5.78 

5.60 

67.7 

>400 

237.1 

273.0 

6 

structure  5  t 
water  ballast  I 

3” 

HM 

HT 

5.58 

5.46 

67.7 

>400 

240.2 

275.0 

7 

structure  5  t 
water  ballast  11 

3" 

HM 

HT 

5.33 

5.24 

67.7 

>400 

238.7 

290.0 

8 

structure  5 

antlsyneetrlc 

modes 

3* 

HM 

HT 

7.0 

67.7 

Table  6  FASTOP  CALCULATION  RESULTS 


FOR  SB13 


Fig.  27  MODE  SHAPES  AND  FLUTTER 

CALCULATION  RESULTS  FOR  THE 
nNAL  FASTOP  DESIGN 


Fig.  29  FLUTTER  CALCULATION  RESULT 
FOR  WATER  BALLAST  CONFIGU¬ 
RATION  II 
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There  was  also  a  shortcoming  in  the  US-swept  forward 
wing  aircraft  X29  which  was  solely  optimized  to  achieve  a 
high  divergence  speed.  With  a  combination  of  aeroelastic  tai¬ 
loring  and  active  control  [17]  high  structural  weight  savings 
can  be  achieved  if  such  a  problem  exists.  Since  CCV  is  quite 
common  now  for  fighter  airplanes  flight  certificarion  for  such 
a  flight  control  system  could  be  received.  Most  important  is 
that  the  lay  out  is  done  in  the  design  stage  and  not  as  a  repair 
solutions. 

It  should  also  be  mentioned  here  that  the  interaction  bet¬ 
ween  rigid  body  and  elastic  structure  shows  the  necessity  to 
inccxporate  flight  mechanics  in  modem  design-  and  optimiza¬ 
tion  programs.  A  twin-engine  prototype  from  Partenavia  was 
lost  in  a  fatal  accident  because  of  the  coupling  between  a  ho¬ 
rizontal  tail  tab-mode  and  the  short  period  motion  [18]. 

The  use  of  aeroelastic  tailoring  makes  it  even  more  im¬ 
portant  to  study  aeroelasticity  parallel  to  other  disciplines  (in 
the  design)  and  not  in  series  as  it  was  done  in  the  past. 

Compared  to  the  initial  design,  the  rtwre  than  1 10%  in¬ 
crease  in  flutter  speed  with  a  small  weight  increase  shows  the 
potential  of  new  carbon  fibers  and  the  use  of  aeroelastic  tailo¬ 
ring. 


Since  1982  the  TSO  program  (Aeroelastic  Tailoring  and 
Structural  Optmization  [19])  is  in  use  at  MBB  as  a  prelimina¬ 
ry  design  tool  for  aerodynamic  surfaces.  In  this  program  the 
structure  of  the  surface  is  represented  as  a  continuous  plate 
with  variable  thickness.Design  variables  are  coefficients,  des¬ 
cribing  the  thickness  distribution  of  different  composite  la¬ 
yers,  the  fiber  orientation,  and,  if  necessary,  variable  concen¬ 
trated  masses  for  fluner  optimization.  Due  to  a  wide  range  of 
aeroelastic  constraints  such  as  frequency,  flutter  speed,  defor¬ 
mations,  aeroelastic  effectiveness  and  divergence  speed,  the 
program  is  very  suitable  for  aeroelastic  tailoring. 

In  1986,  TSO  was  used  for  a  design  study  on  a  light 
combat  aircraft  wing  shown  in  Fig.  30.  In  this  case,  three  sta¬ 
tic  load  cases  were  used  for  the  pure  strength  design  with  the 
main  objective:  minimum  weight. 

The  flexible  wing  toll  rate  was  not  a  constraint  in  the  be- 
ginning.The  wing  flap  hinge  moments,  however,  are  often  a 
critical  design  criterium  that  can  also  limit  the  toll  rate  of  the 
aircraft.  It  could  be  demonstrated  in  this  study,  that  by  defi¬ 
ning  an  aileron  toll  moment  effectiveness  constraint,  the  wing 
cover  skin  thickness  and  fiber  orientations  could  be  designed 
for  higher  roll  rates  and  a  considerable  reduction  in  hinge  nx>- 
ment  with  a  small  weight  increase. 

The  relations  between  structural  weight,  roll  effective¬ 
ness  and  flap  hinge  moment  can  be  seen  in  Fig.  31.  With  very 
little  increase  in  weight  the  hinge  moment  for  the  requited  toll 
rate  of  180°  at  maximum  dynamic  pressure  and  Mach  1.1  can 
be  reduced  to  30%  of  the  original  value.  This  sensitivity  ana¬ 
lysis  also  indicates  the  region  where  an  additional  increase  in 
structural  mass  can  not  improve  the  performance.  Without  the 
use  of  an  optimization  program,  this  kind  of  trend  studies 
would  be  impossible,  especially  if  composite  materials  or  ae- 
Toelastic  requirements  ate  involved.  Fig.  32  gives  a  compari¬ 
son  of  the  wing  cover  skin  designs  for  strength  requirements 
only  and  for  an  additional  flexible  aileron  roll  effectiveness 
requirementThe  importance  of  the  seletion  of  design  varia¬ 


bles  is  demonstrated  in  Fig.  33  for  four  different  approaches: 


fixed  fiber  orientation  and 
balanced  -45°  and  -45°  plies 

W 

=  100% 

fixed  fiber  orientation, 
unbalanced 

w 

=  76% 

free  fiber  orientation, 
balanced 

w 

=  67% 

flee  fiber  orientation, 
unbalaiKed 

w 

=  54% 

Fig,  30  LIGHT  COMBAT  AIRCRAFT  WING 
PLATE  MODEL 


Fig.  31  OPTIMIZATION  WING  COVER 
SKIN  WEIGHT  FOR  DIFFERENT 
ROLL  EFFECTIVENESS  CON¬ 
STRAINTS  AND  REQUIRED  FLAP 
HINGE  MOMENT  FOR  180  DEG/SEC 
ROLL  RATE  AT  MA.  1.1,  SEA  LEVEL 


Hinge  Moment  kNm 
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Fig.  32  OPTIMIZED  WING  SKIN  THICK¬ 
NESS  DISRIBUTIONS  FOR  DESIGNS 
WITH  AND  WITHOUT  FLAP  ROLL 
MOMENT  EFFECTIVENESS  CON¬ 
STRAINT 


Fig.  33  FOUR  DIFFERENT  DESIGNS  FOR 
STRENGTH  AND  ROLL  EFFECTI¬ 
VENESS  CONSTRAINT  (RM). 
ARROWS  INDICATE  AMOUNT  OF 
PLIES  IN  EACH  DIRECTION 


APYANCEB-ABCKAP  PEStGN  IS3TH 

MBB-LAGRANGE 

A  large  number  of  other  studies  and  applications  of  LA¬ 
GRANGE  to  current  projects  have  been  performed  already 
and  presented  in  several  publications  [20],  [21],  [22], 

A  typical  example  for  the  application  of  LAGRANGE  is 
the  composite  wing  structure  for  the  experimental  aircraft  X- 
31  A.  In  this  case  optimization  was  beneficial  for  two  main 
objectives  of  the  program;  a  low  cost  approach  and  a  very 
sh<xt  tune  for  develtqiment  and  design.  Besides  a  design  for 
minimum  weight  another  requirement  was  a  high  flutter  mar¬ 
gin  to  reduce  efforts  and  costs  for  flutter  wind  tunnel  and 
flight  tests.  Although  flutter  did  not  effea  the  design,  it  could 
be  surveyed  simultaneously  during  c^timization.  Static  aero- 
elastic  effectiveness  was  also  investigated  during  the  design 
process. 


Fig.  34  FINITE  ELEMENT  MODEL  OF  THE 
X-31  WING 


A  finite  element  model  of  the  wing  is  depicted  in  Fig.  34. 
It  has  1764  elements  and  1871  degrees  of  fleedom.  The  opti¬ 
mized  sldn  thicknesses  were  then  translated  into  design  dra¬ 
wings  with  small  modifications.  As  an  example,  the  upper 
wing  sldn  weight  of  S3  kg  flora  an  initial  design  (preoptimized 
with  another  program)  could  be  reduced  to  44  kg  in  the  FEM 
which  resulted  in  4S  kg  in  the  actual  design.  The  final  design 
meets  the  target  weight  and  has  a  margin  of  100%  in  airspeed 
for  flutter. 

INTEGRATED  DESIGN  CONCEPT 

The  experience  obtained  from  various  designs  with  design 
requirements  coming  from  different  disciplines  has  shown  the 
ne^  for  integrated  design  concepts  and  programs. 

The  interactions  between  aerodynamics,  performance, 
flight  control,  structural  loads,  dynamics,  aeroelasticity, 
strength  and  materials,  and  finally  the  design  have  always  exi¬ 
sted. 

Due  to  increasing  aitcrafl  performance  requirements,  like 
payload,  fuel  efficiency,  or  maneuvering  capability,  these  inte¬ 
ractions  have  become  stronger  and  more  imponant 

The  mass  is  very  important  for  an  airciaft  The  aerodyna¬ 
mic  lift  or  drag  from  different  sotirces  like  surfoce  area,  distri¬ 
bution  of  cross  section  areas,  lift-induced  drag  and  airfbil  sha¬ 
pe  ate  also  main  design  parameters.  The  influence  of 
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geometric  parameters  of  an  aircraft  on  lift  and  drag  characte¬ 
ristics  has  been  studied  since  the  earliest  days  of  aviation. 

The  importance  of  aspect  ratio  on  lift  induced  drag,  of 
wing  taper  ratio  on  lift  distribution  and  drag  of  wing  taper  ra¬ 
tio  on  lift  distribution  and  drag,  of  wing  thickness  ratio  and 
sweep  angle  on  drag  increase  with  Mach  number  is  well 
known.  But  how  is  the  influence  of  these  parameters  on  the 
structural  mass  and  on  the  loads  that  cause  the  mass? 


To  demonstrate  these  relations,  parametric  studies  have 
been  performed  for  typical  fighter  aircraft  wings,  using  the 
TSO  proram.  Some  of  these  wings  are  depicted  in  Fig.  35. 

of. 


S«««p  Anqlt 


Aspect  lUtle 


T«p«r  Ast^o 


Fig.  35  WING  GEOMETRY  PARAMETERS 
FOR  OPTIMIZATION  STUDIES 

For  all  these  wings  the  same  basic  design  requirements 
have  been  defined,  using  identical  total  weight  and  the  same 
wing  area: 

•  a  9  g  static  load  case  with  aetoelastically  trimmed  load 
distribution 

•  a  maximum  roll  rate  at  high  dynamic  pressure  for  the  ai¬ 
leron  effectiveness 

•  a  minimum  flutter  speed  of  1000  kts 


The  carbon  fiber  wing  cover  skins  of  these  wings  were 
then  optimized  for  the  above  load  cases  separately  and  simul¬ 
taneously. 

Fig.  36  depicts  the  influence  of  the  aspect  ratio  on  the 
skin  weight  for  the  different  design  requirements.  Of  course, 
a  change  in  requirements,  a  different  flap  geometry  or  a  diffe¬ 
rent  mass  distribution  will  give  different  trends  for  the  geo¬ 
metry  parameters  and  for  the  additional  weight  required  to 
meet  each  individual  requirement  on  lop  of  others. 

Therefore,  similar  studies  should  be  performed  in  the  pre¬ 
liminary  stage  of  a  new  aircraft  design. 


Fig.  36  OPTIMIZED  WING  COVER  SKIN 
WEIGHT  VS  ASPECT  RATIO 


AEROELASTIC  TAILORING  OF  A  FIN  MADE  OF 
COMPOSITE  MATERIAL  1231 

An  aircraft  fin  has  to  fullfil  quite  different  design  require¬ 
ments  with  a  similar  priority  and  the  flnal  design  requires  the 
evaluation  of  many  off-design  point  studies. 

The  design  of  aerodynamic  surfaces  such  as  wing,  fin,  fo¬ 
replane  and  tailplane  needs  two  major  design  steps: 

First,  the  aerodynamic  design  to  define  the  overall  geo¬ 
metry  like  area,  span,  aspect  ratio,  uper  ratio  and  profit. 

Second,  the  structural  design  to  develop  the  internal  scmc- 
tural  arrangement  of  skin,  ribs,  stringers,  spars,  rudder  sup¬ 
port,  rudder  acmation,  anachments,  equipment  systems. 

The  final  design  must  fulfill  the  following  design  require¬ 
ments  with  minimum  weight: 

•  Static  strength  to  withstand  design  loads 

•  Aeroelastic  efficiencies  for  performance 

•  No  flutter  inside  of  the  flight  envelope 


Fig.  37  FIN  STRUCTURAL  MODEL  WITH 
SKIN  ELEMENT  NUMBERS 
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Ii  is  quite  clear  that  such  a  design  requires  an  interactive 
coupling  of  the  above  mentioned  two  design  steps.  A  structu¬ 
ral  model  of  the  investigated  fin  is  shown  in  Fig.  37..  A  com¬ 
parison  of  initial  design  analysis  results  and  design  constraints 
is  given  in  the  following  table; 


strength 

DESIGN 

CONSTRAINT 

St ( a  1 n  all owab t « & 

Tent  1  on  <  . 00  3  T 
Conpariton  5-. 0028 

INITIAL 

DESIGN 

Loa(3  rase  2 

E 1 eaen  t  1 6 

Eleaent  22 

VIOLATION 

-  123 

-  228 

AEROELASTICS 

efficiencies 
riN  .8 

.  TS3 

-  .  059 

'..8  'bOFTS 

RUDDER  .S 

.441 

•  .  118 

FLvTTEP 

Ma  : . 2/S . L. 

Flutter  speed 

vr  -  530  B  sec 

495  pi  sec 

-  .  066 

The  frequency  versus  speed  behaviour  for  the  optimi¬ 
zed/initial  structure  is  given  in  Fig.  38.  The  corresponding 
damping  is  plotted  in  Fig.  39.  The  results  of  the  optimization 
procedure  are  shown  in  Table  8. 

The  flutter  speed  is  increased  to  530  m/sec.  anu  aeroela- 
stic  efficiencies  are  increased  8%  for  the  fin  and  for  the  rud¬ 
der  by  13%.  The  structural  weight  is  reduced  by  7%. 

Skin  thicknesses  for  the  different  carbon  fibre  layer  of 
the  optimum  strucmral  design  are  presented  in  Fig.  40-43. 


MBB-LRGRONGE 

FLUTTER  FREQUENCY  PLOT 
RCB-FIN  tOPTlMIZEOl 

sTnTic/pepaeLNtttc/FLUTTeM 


Optimiaed  Design:  _  V,,_  -  530  m/s 


Fig.  38  FLUTTER  ANALYSIS  FREQUENCY 
PLOT  FOR  INITIAL  AND  OPTIMI¬ 
ZED  DESIGN 


MBB-LPIGRflNGE 

FLUTTER  ORMPINO  PLOT 
RCR-FIN  (OPTIMIZED) 


Fig.  39  FLUTTER  ANALYSIS  DAMPING 
PLOT  FOR  INITIAL  AND  OPTIMI¬ 
ZED  DESIGN 


DESIGN 

INITIAL 

CCNSTRAINT 

VIOLATED 

DESICN 

OPTIMAL 

WEIGHT 

[kg] 

Structure 

99.4 

92.9 

Non  Struc¬ 
ture 

53.6 

53.6 

Total 

153. 

146.5 

STREICTH 

Loadcase  2 
Element 

18  -  .123 

22  -  .224 

-.123 

-.228 

EIGEN- 

FREQUorcy 

[Hz] 

fl  -  8.90 
f2  -  29.83 
f3  -  31.16* 
f4  -  39.97 
f5  -  54.86 

9.20 

30.21(x-l.) 

30.61 

41.08 

58.31 

FIAOTER 

FREQUENCY 

fF  .  21.22 

22.0 

SPEED  Im/S) 

VF  -  495. 

-.066 

530. 

AERO¬ 

ELASTICS 

FIN 

RODDER 

.753 

.441 

-.059 

-.118 

.814 

.500 

*  (x-l.) 

Table  8  OPTIMIZATION  RESULTS 
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Fig.  40  AEROELASTIC  DESIGN:  SKIN 
THICKNESS  FOR  LAYER  1 


Fig.  41  AEROELASTIC  DESIGN:  SKIN 
THICKNESS  FOR  LAYER  2 


Fig.  42  AEROELASTIC  DESIGN :  SKIN 
THICKNESS  FOR  LAYER  3 


Fig.  43  AEROELASTIC  DESIGN:  SKIN 
THICKNESS  FOR  LAYER  4 


CURRENT  ACTIVITIES  AT  MBB  MUNICH  IN 
THE  FIELD  OF  AEROELASTIC  TAILORING 

Shane  Optimization: 

For  aeroelastic  applications,  the  variation  of  the  external 
geometry  is  most  important.  The  problem  in  dealing  with  grid 
point  coordinates  as  design  variables  is  the  connection  with 
the  aerodynamic  model. 

The  variation  of  the  aerodynamic  model  and  its  elements 
has  to  be  combined  with  the  structural  model.  This  synthesis 
is  currently  being  investigated. 

Topology  Optimiiatioii 

As  a  higher  level  aspect  of  structural  optimization  pro¬ 
blems  (after  sizing  and  shape  optimization),  topology  optimi¬ 
zation  of  the  internal  and  external  structural  geometry  is  pio- 
missing  great  potentials  in  aircraft  perfwmance  and  structural 
design.  Based  on  experience  gained  during  a  study  which  has 
been  made  to  find  the  optimal  attachment  coordinates  for 
wing-flap  connections  and  with  basic  analytical  tools  develo¬ 
ped  for  general  investigations  of  topological  aspects,  an  ex¬ 
tended  version  of  this  program  is  currently  planned  for  air¬ 
craft  sBuctures. 

Smart  Acttwtors 

A  new  type  of  actuators  with  multi-signal  input  capability 
is  currently  under  development.  This  actuator  will  offer  seve¬ 
ral  advantages  in  aeroelasticity:  it  will  improve  the  aeioser- 
voelastic  stability  of  the  system,  it  could  be  used  to  replace 
buzz  dampers,  and  it  is  capable  to  cover  active  control  tech¬ 
nology  for  aeroelastic  aspects  without  additional  wotlc  load 
for  the  flight  control  computer. 

AIRCRAFT  OPTIMIZATION  PROGRAM 

Based  on  MBB-LAGRANGE,  an  optimization  program 
for  a  total  aircraft  is  currently  being  developed.  In  a  first  step, 
it  will  integrate  the  external  loads  in  the  optimization  process. 
Because  the  complete  aircraft  structure  is  used,  trimmed  ftee- 
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free  conditions  will  be  simulated  to  describe  the  coirect  load 
distribution,  inlcuding  static  aeroelastic  defcHmations.  Other 
aeroelastic  problems  like  flutter  and  gust  response  will  also 
be  analysed  and  included  in  the  t^timization  process  for  the 
complete  conriguration  from  the  beginning  of  the  design  pro¬ 
cess. 

The  program  will  then  be  extended  to  cover  aerodynamic 
aspects  like  drag,  UD,  Cj^,  optimum  twist  and  camber  or 
design  lift  conditions. 

For  this  purpose,  the  aerodynamic  model  of  the  aircraft 
will  be  used  to  provide  design  variables  in  addition  to  the 
structural  variables.  This  combination  will  also  allow  to  cover 
other  important  aspects  like  stability,  performance  and  con¬ 
trol  and  it  will  provide  a  basis  for  the  integration  of  active 
control  techniques  into  the  optimization  loop. 

The  challenge  on  aeroelastic  tailoring  is  depending  on 
new  flight  performance  requirements,  demanding  new  confi¬ 
gurations  as  well  as  new  technologies  and  new  materials. 

At  the  present  time  the  development  of  new  airplanes  is 
influenced  by  new  techniques,  such  as  flutter  suppression, 
CCV-configuration,  gust  load  alleviation  etc.  In  addition  to 
stress,  displacement,  aeroelastic  and  dynamic  constraints  an 
integrated  design  involves  all  these  techniques  and  the  opti¬ 
mization  procedures  must  be  extended  for  these  new  con¬ 
straints. 

Especially  the  combination  of  new  developments  in  aero¬ 
dynamic  shape  optimization  and  the  well  experienced  active 
control  technologies  with  structural  optimization  routines  will 
necessarily  enter  into  a  multidisciplinary  optimization  pro¬ 
cess. 

It  will  take  a  further  period  of  developtnent  even  if  the 
progress  in  computer  power  and  new  mathematical  optimiza¬ 
tion  methods  are  enormous.  Existing  technologies  have  to  be 
refined,  new  developments  like  shape  optimization  still  have 
to  be  completed  and  experiences  gained  with  this  tools,  and, 
last  not  least,  specialists  of  different  disciplines  have  to  be 
convinced  that  the  new  opportunities  are  worth  the  effort  ta¬ 
king  into  the  bargain  a  highly  increased  coitqrlexity  of  the  de¬ 
sign  process.  To  account  for  the  increased  corttpiexity,  adap¬ 
ted  intelligent  user  interfaces  and  checking  routines  for  the 
generation  of  reliable  inputs,  for  the  check  of  interim  solu¬ 
tions  and  for  the  interpretation  of  the  output  as  well  as  impro¬ 
ved  integrated  expert  systems  to  suppmt  the  selection  of  ap- 
propiaie  algorithm  are  required. 

Up  to  this  status  it  can  not  be  expected  that  there  will  be 
one  program  only  for  the  optimal  design  during  different  sta¬ 
ges  of  aircraft  projects. 

An  initial  preliminary  nrsign  should  in  fact  include  as 
many  desciplines  as  possibles.  But  at  the  same  time,  this  task 
must  remain  in  a  not  too  detailed  and  complex  level  to  allow 
the  investigation  of  a  great  number  of  designs  and  to  answer 
questions  concerning  essential  changes  of  design  require¬ 
ments  in  a  relatively  short  period  of  time.  After  this,  the  indi¬ 
vidual  disciplines  should  use  their  own  programs  and  me¬ 
thods  to  find  the  tqytimum  in  a  more  detailed  model,  without 
forgetting  the  neighbour  areas. 

The  preliminary  design  program  could  in  parallel  serve  as 
a  tool  to  integrate  the  results  from  detail  designs. 

Large  efforts  will  be  required  to  reduce  the  enoimous 
computational  costs  by  the  development  of  efficient  methods 
for  cross  sensitivity  calculations  and  for  approximate  optimi¬ 


zation  procedures. 

An  extensive  descriptibn  of  the  use  of  cqitimization  for 
concurrent  engineering  is  given  in  [24]. 
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1 4.  Abstract 

There  is  a  lack  of  precise  information  on  the  effectiveness  of  specific  methods  in  generating 
optimum  designs  for  realistic  aircraft  structures.  In  this  situation  it  is  difficult  for  designers  to 
make  decisions  on  which  systems  to  employ  for  a  given  design  problem  and  which  developments 
to  pursue.  Thus  it  is  neces.sary  for  designers  to  be  aware  of  the  relative  merits  of  the  different 
methods  currently  used  for  the  design  optimisation  of  advanced  aircraft. 

This  Lecture  Series  covers  the  methods  available  for  the  computer  based  design  analysis  and 
design  optimisation  of  aircraft  structures.  The  Lecture  Series  deals  with  the  principles  and 
practices  adopted  to  integrate  the  various  factors  which  are  considered  in  the  design  of  advanced 
aircraft.  These  factors  include:  structural  shape,  aerodynamics,  active  control  technology  and 
aircraft  performance.  Realistic  case  studies  arc  used  to  illustrate  the  methods  used  for  different 
design  problems. 

The  following  topics  are  covered  in  detail: 

—  Overview  of  integrated  design  analysis,  background,  methods,  objectives  and  requirements. 

—  Optimisation  in  design  (CAE/CAD). 

—  A  system  approach  to  aircraft  optimisation. 

—  Case  studies  for  different  design  proulems. 


This  Lecture  Series,  sponsored  by  the  Structures  and  Materials  Panel  of  AGARD.  has  been 
implemented  by  the  Consultant  and  Exchange  Programme. 
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